1package Graph;
2
3use strict;
4use warnings;
5BEGIN { warnings->unimport('recursion') if $ENV{GRAPH_ALLOW_RECURSION} }
6
7sub __carp_confess { require Carp; Carp::confess(@_) }
8BEGIN {
9    if (0) { # SET THIS TO ZERO FOR TESTING AND RELEASES!
10	$SIG{__DIE__ } = \&__carp_confess;
11	$SIG{__WARN__} = \&__carp_confess;
12    }
13}
14
15use Graph::AdjacencyMap qw(:flags :fields);
16
17our $VERSION = '0.9725';
18
19require 5.006; # Weak references are absolutely required.
20
21my @GRAPH_PROPS_COPIED = qw(
22    undirected refvertexed countvertexed multivertexed __stringified
23    hyperedged countedged multiedged
24);
25my $_empty_array = [];
26sub _empty_array () { $_empty_array }
27
28my $can_deep_copy_Storable;
29sub _can_deep_copy_Storable () {
30    return $can_deep_copy_Storable if defined $can_deep_copy_Storable;
31    return $can_deep_copy_Storable = 0 if $] < 5.010; # no :load tag Safe 5.8
32    eval {
33        require Storable;
34        require B::Deparse;
35        Storable->VERSION(2.05);
36        B::Deparse->VERSION(0.61);
37    };
38    $can_deep_copy_Storable = !$@;
39}
40
41sub _F () { 0 } # Flags.
42sub _G () { 1 } # Generation.
43sub _V () { 2 } # Vertices.
44sub _E () { 3 } # Edges.
45sub _A () { 4 } # Attributes.
46sub _U () { 5 } # Union-Find.
47
48my $Inf;
49
50BEGIN {
51  if ($] >= 5.022) {
52    $Inf = eval '+"Inf"'; # uncoverable statement
53  } else {
54    local $SIG{FPE}; # uncoverable statement
55    eval { $Inf = exp(999) } || # uncoverable statement
56	eval { $Inf = 9**9**9 } || # uncoverable statement
57	    eval { $Inf = 1e+999 } || # uncoverable statement
58		{ $Inf = 1e+99 }; # uncoverable statement
59                # Close enough for most practical purposes.
60  }
61}
62
63sub Infinity () { $Inf }
64
65# Graphs are blessed array references.
66# - The first element contains the flags.
67# - The second element is the vertices.
68# - The third element is the edges.
69# - The fourth element is the attributes of the whole graph.
70# The defined flags for Graph are:
71# - unionfind
72# The vertices are contained in a "simplemap"
73# (if no attributes) or in a "map".
74# The edges are always in a "map".
75# The defined flags for maps are:
76# - _COUNT for countedness: more than one instance
77#   expects one for vertices and two for edges
78# - _UNORD for unordered coordinates (a set): if _UNORD is not set
79#   the coordinates are assumed to be meaningfully ordered
80# Vertices and edges assume none of these flags.
81
82use Graph::Attribute array => _A, map => 'graph';
83
84sub stringify {
85    my ($u, $h) = (&is_undirected, &is_hyperedged);
86    my $e = $u ? '=' : '-';
87    my @edges = map join($e,
88	$u ? sort { "$a" cmp "$b" } @$_ :
89	$h ? map '['.join(",", sort { "$a" cmp "$b" } @$_).']', @$_ :
90	@$_), &_edges05;
91    my @s = sort @edges;
92    push @s, sort { "$a" cmp "$b" } &isolated_vertices;
93    join(",", @s);
94}
95
96sub eq {
97    "$_[0]" eq "$_[1]"
98}
99
100sub boolify {
101  1;  # Important for empty graphs: they stringify to "", which is false.
102}
103
104sub ne {
105    "$_[0]" ne "$_[1]"
106}
107
108use overload
109    '""' => \&stringify,
110    'bool' => \&boolify,
111    'eq' => \&eq,
112    'ne' => \&ne;
113
114sub _opt {
115    my ($opt, $flags, %flags) = @_;
116    while (my ($flag, $FLAG) = each %flags) {
117	$$flags |= $FLAG if delete $opt->{$flag};
118	$$flags &= ~$FLAG if delete $opt->{"non$flag"};
119    }
120}
121
122sub _opt_get {
123    my ($opt, $key, $var) = @_;
124    return if !exists $opt->{$key};
125    $$var = delete $opt->{$key};
126}
127
128sub _opt_unknown {
129    my ($opt) = @_;
130    return unless my @opt = keys %$opt;
131    __carp_confess sprintf
132        "@{[(caller(1))[3]]}: Unknown option%s: @{[map qq['$_'], sort @opt]}",
133        @opt > 1 ? 's' : '';
134}
135
136sub _opt_from_existing {
137    my ($g) = @_;
138    my %existing;
139    $existing{$_}++ for grep $g->$_, @GRAPH_PROPS_COPIED;
140    $existing{unionfind}++ if $g->has_union_find;
141    %existing;
142}
143
144sub _opt_to_vflags {
145    my ($vflags, $opt) = (0, @_);
146    _opt($opt, \$vflags,
147	 countvertexed	=> _COUNT,
148	 multivertexed	=> _MULTI,
149	 refvertexed	=> _REF,
150	 refvertexed_stringified => _REFSTR ,
151	 __stringified => _STR,
152	);
153    $vflags;
154}
155
156sub _opt_to_eflags {
157    my ($eflags, $opt) = (0, @_);
158    $opt->{undirected} = !delete $opt->{directed} if exists $opt->{directed};
159    _opt($opt, \$eflags,
160	 countedged	=> _COUNT,
161	 multiedged	=> _MULTI,
162	 undirected	=> _UNORD,
163	);
164    ($eflags, delete $opt->{hyperedged});
165}
166
167sub new {
168    my ($class, @args) = @_;
169    my $gflags = 0;
170    my %opt = _get_options( \@args );
171
172    %opt = (_opt_from_existing($class), %opt) # allow overrides
173	if ref $class && $class->isa('Graph');
174
175    my $vflags = _opt_to_vflags(\%opt);
176    my ($eflags, $is_hyper) = _opt_to_eflags(\%opt);
177
178    _opt(\%opt, \$gflags,
179	 unionfind     => _UNIONFIND,
180	);
181
182    my @V;
183    if ($opt{vertices}) {
184	__carp_confess "Graph: vertices should be an array ref"
185	    if ref $opt{vertices} ne 'ARRAY';
186	@V = @{ delete $opt{vertices} };
187    }
188
189    my @E;
190    if ($opt{edges}) {
191        __carp_confess "Graph: edges should be an array ref of array refs"
192	    if ref $opt{edges} ne 'ARRAY';
193	@E = @{ delete $opt{edges} };
194    }
195
196    _opt_unknown(\%opt);
197
198    __carp_confess "Graph: both countvertexed and multivertexed"
199	if ($vflags & _COUNT) && ($vflags & _MULTI);
200
201    __carp_confess "Graph: both countedged and multiedged"
202	if ($eflags & _COUNT) && ($eflags & _MULTI);
203
204    my $g = bless [ ], ref $class || $class;
205
206    $g->[ _F ] = $gflags;
207    $g->[ _G ] = 0;
208    $g->[ _V ] = _make_v($vflags);
209    $g->[ _E ] = _make_e($is_hyper, $eflags);
210
211    $g->add_vertices(@V) if @V;
212
213    __carp_confess "Graph: edges should be array refs"
214	if grep ref $_ ne 'ARRAY', @E;
215    $g->add_edges(@E);
216
217    $g->[ _U ] = do { require Graph::UnionFind; Graph::UnionFind->new }
218	if $gflags & _UNIONFIND;
219
220    return $g;
221}
222
223sub _make_v {
224    my ($vflags) = @_;
225    $vflags ? _am_heavy($vflags, 1) : _am_light($vflags, 1);
226}
227
228sub _make_e {
229    my ($is_hyper, $eflags) = @_;
230    ($is_hyper or $eflags & ~_UNORD) ?
231	_am_heavy($eflags, $is_hyper ? 0 : 2) :
232	    _am_light($eflags, 2);
233}
234
235sub _am_light {
236    require Graph::AdjacencyMap::Light;
237    Graph::AdjacencyMap::Light->_new(@_);
238}
239
240sub _am_heavy {
241    Graph::AdjacencyMap->_new(@_);
242}
243
244sub countvertexed { $_[0]->[ _V ]->_is_COUNT }
245sub multivertexed { $_[0]->[ _V ]->_is_MULTI }
246sub refvertexed   { $_[0]->[ _V ]->_is_REF   }
247sub refvertexed_stringified { $_[0]->[ _V ]->_is_REFSTR }
248sub __stringified { $_[0]->[ _V ]->_is_STR   }
249
250sub countedged    { $_[0]->[ _E ]->_is_COUNT }
251sub multiedged    { $_[0]->[ _E ]->_is_MULTI }
252sub hyperedged    { !$_[0]->[ _E ]->[ _arity ] }
253sub undirected    { $_[0]->[ _E ]->_is_UNORD }
254
255sub directed { ! $_[0]->[ _E ]->_is_UNORD }
256
257*is_directed      = \&directed;
258*is_undirected    = \&undirected;
259
260*is_countvertexed = \&countvertexed;
261*is_multivertexed = \&multivertexed;
262*is_refvertexed   = \&refvertexed;
263*is_refvertexed_stringified = \&refvertexed_stringified;
264
265*is_countedged    = \&countedged;
266*is_multiedged    = \&multiedged;
267*is_hyperedged    = \&hyperedged;
268
269sub has_union_find { $_[0]->[ _U ] }
270
271sub add_vertex {
272    __carp_confess "Graph::add_vertex: use add_vertices for more than one vertex" if @_ != 2;
273    __carp_confess "Graph::add_vertex: undef vertex" if grep !defined, @_;
274    goto &add_vertices;
275}
276
277sub has_vertex {
278    my $g = $_[0];
279    my $V = $g->[ _V ];
280    return defined $V->has_path($_[1]) if ($V->[ _f ] & _REF);
281    exists $V->[ _pi ]->{ $_[1] };
282}
283
284sub _vertices05 {
285    my $g = $_[0];
286    $g->[ _V ]->paths;
287}
288
289sub vertices {
290    my $g = $_[0];
291    my @v = &_vertices05;
292    return @v if !(&is_multivertexed || &is_countvertexed);
293    return map +(($_) x $g->get_vertex_count($_)), @v if wantarray;
294    my $V = 0;
295    $V += $g->get_vertex_count($_) for @v;
296    return $V;
297}
298
299*unique_vertices = \&_vertices05;
300
301sub has_vertices {
302    my $g = shift;
303    scalar $g->[ _V ]->has_any_paths;
304}
305
306sub add_edge {
307    &expect_hyperedged, &expect_undirected if @_ != 3;
308    $_[0]->add_edges([ @_[1..$#_] ]);
309}
310
311sub _vertex_ids_ensure {
312    push @_, 1;
313    goto &_vertex_ids_maybe_ensure;
314}
315
316sub _vertex_ids_ensure_multi {
317    my $id = pop;
318    my @i = &_vertex_ids_ensure;
319    push @_, $id;
320    @i ? (@i, $id) : ();
321}
322
323sub _vertex_ids {
324    push @_, 0;
325    goto &_vertex_ids_maybe_ensure;
326}
327
328sub _vertex_ids_multi {
329    my $id = pop;
330    my @i = &_vertex_ids;
331    push @_, $id;
332    @i ? (@i, $id) : ();
333}
334
335sub _vertex_ids_maybe_ensure {
336    my $ensure = pop;
337    my ($g, @args) = @_;
338    __carp_confess "Graph: given undefined vertex" if grep !defined, @args;
339    my $V = $g->[ _V ];
340    my $deep = &is_hyperedged && &is_directed;
341    return $V->get_ids_by_paths(\@args, $ensure, $deep) if ($V->[ _f ] & _REF) or $deep;
342    my $pi = $V->[ _pi ];
343    my @non_exist = grep !exists $pi->{ $_ }, @args;
344    return if !$ensure and @non_exist;
345    $V->get_ids_by_paths(\@non_exist, 1) if @non_exist;
346    @$pi{ @args };
347}
348
349sub has_edge {
350    my $g = $_[0];
351    my $E = $g->[ _E ];
352    my ($Ef, $Ea) = @$E[ _f, _arity ];
353    return 0 if $Ea and @_ != $Ea + 1;
354    my $directed = &is_directed;
355    my $deep = &is_hyperedged && $directed;
356    return 0 if (my @i = &_vertex_ids) != @_ - 1;
357    return defined $E->has_path($directed ? \@i : [ map [ sort @$_ ], @i ]) if $deep;
358    @i = sort @i if !$directed;
359    exists $E->[ _pi ]{ "@i" };
360}
361
362sub any_edge {
363    my ($g, @args) = @_;
364    my $E = $g->[ _E ];
365    my $V = $g->[ _V ];
366    return 0 if (my @i = $V->get_ids_by_paths(\@args)) != @args;
367    $E->has_successor(@i);
368}
369
370sub _edges05 {
371    my $g = $_[0];
372    my @e = $g->[ _E ]->paths;
373    return @e if !wantarray;
374    $g->[ _V ]->get_paths_by_ids(\@e, &is_hyperedged && &is_directed);
375}
376
377*unique_edges = \&_edges05;
378
379sub edges {
380    my $g = $_[0];
381    my @e = &_edges05;
382    return @e if !(&is_multiedged || &is_countedged);
383    return map +(($_) x $g->get_edge_count(@$_)), @e if wantarray;
384    my $E = 0;
385    $E += $g->get_edge_count(@$_) for @e;
386    return $E;
387}
388
389sub has_edges {
390    scalar $_[0]->[ _E ]->has_any_paths;
391}
392
393###
394# by_id
395#
396
397sub add_vertex_by_id {
398    &expect_multivertexed;
399    my ($g, $v, $id) = @_;
400    my $V = $g->[ _V ];
401    return $g if $V->has_path_by_multi_id( my @args = ($v, $id) );
402    my ($i) = $V->set_path_by_multi_id( @args );
403    $g->[ _U ]->add($i) if &has_union_find;
404    $g->[ _G ]++;
405    return $g;
406}
407
408sub add_vertex_get_id {
409    &expect_multivertexed;
410    my ($g, $v) = @_;
411    my ($i, $multi_id) = $g->[ _V ]->set_path_by_multi_id( $v, _GEN_ID );
412    $g->[ _U ]->add($i) if &has_union_find;
413    $g->[ _G ]++;
414    return $multi_id;
415}
416
417sub has_vertex_by_id {
418    &expect_multivertexed;
419    my ($g, $v, $id) = @_;
420    $g->[ _V ]->has_path_by_multi_id( $v, $id );
421}
422
423sub delete_vertex_by_id {
424    &expect_multivertexed;
425    &expect_non_unionfind;
426    my ($g, $v, $id) = @_;
427    return $g unless &has_vertex_by_id;
428    # TODO: what to about the edges at this vertex?
429    # If the multiness of this vertex goes to zero, delete the edges?
430    $g->[ _V ]->del_path_by_multi_id( $v, $id );
431    $g->[ _G ]++;
432    return $g;
433}
434
435sub get_multivertex_ids {
436    &expect_multivertexed;
437    my $g = shift;
438    $g->[ _V ]->get_multi_ids( @_ );
439}
440
441sub add_edge_by_id {
442    &expect_multiedged;
443    my $g = $_[0];
444    my @i = &_vertex_ids_ensure_multi;
445    my $id = pop @i;
446    @i = sort @i if &is_undirected;
447    $g->[ _E ]->set_path_by_multi_id( \@i, $id );
448    $g->[ _G ]++;
449    $g->[ _U ]->union(\@i) if &has_union_find;
450    return $g;
451}
452
453sub add_edge_get_id {
454    &expect_multiedged;
455    my $g = $_[0];
456    my @i = &_vertex_ids_ensure;
457    @i = sort @i if &is_undirected;
458    my (undef, $id) = $g->[ _E ]->set_path_by_multi_id( \@i, _GEN_ID );
459    $g->[ _G ]++;
460    $g->[ _U ]->union(\@i) if &has_union_find;
461    return $id;
462}
463
464sub has_edge_by_id {
465    &expect_multiedged;
466    my $g = $_[0];
467    my @i = &_vertex_ids_multi;
468    return 0 if @i < @_ - 2;
469    my $id = pop @i;
470    @i = sort @i if &is_undirected;
471    $g->[ _E ]->has_path_by_multi_id( \@i, $id );
472}
473
474sub delete_edge_by_id {
475    &expect_multiedged;
476    &expect_non_unionfind;
477    my $g = $_[0];
478    my $E = $g->[ _E ];
479    my @i = &_vertex_ids_multi;
480    return if @i < @_ - 2;
481    my $id = pop @i;
482    @i = sort @i if &is_undirected;
483    return unless $E->has_path_by_multi_id( my @args = (\@i, $id) );
484    $E->del_path_by_multi_id( @args );
485    $g->[ _G ]++;
486    return $g;
487}
488
489sub get_multiedge_ids {
490    &expect_multiedged;
491    return unless @_-1 == (my @i = &_vertex_ids);
492    $_[0]->[ _E ]->get_multi_ids( \@i );
493}
494
495###
496# Neighbourhood.
497#
498
499sub _edges_at {
500    goto &_edges_from if &is_undirected;
501    require Set::Object;
502    Set::Object->new(&_edges_from, &_edges_to)->${ wantarray ? \'members' : \'size' };
503}
504
505sub _edges_from {
506    my ($g, @args) = @_;
507    my ($V, $E) = @$g[ _V, _E ];
508    return if (my @i = $V->get_ids_by_paths(\@args, &is_hyperedged && &is_directed)) != @args;
509    $E->paths_from(@i);
510}
511
512sub _edges_to {
513    goto &_edges_from if &is_undirected;
514    my ($g, @args) = @_;
515    my ($V, $E) = @$g[ _V, _E ];
516    return if (my @i = $V->get_ids_by_paths(\@args, &is_hyperedged && &is_directed)) != @args;
517    $E->paths_to(@i);
518}
519
520sub edges_at {
521    goto &_edges_at if !wantarray;
522    $_[0]->[ _V ]->get_paths_by_ids([ &_edges_at ], &is_hyperedged && &is_directed);
523}
524
525sub edges_from {
526    goto &_edges_from if !wantarray;
527    $_[0]->[ _V ]->get_paths_by_ids([ &_edges_from ], &is_hyperedged && &is_directed);
528}
529
530sub edges_to {
531    goto &edges_from if &is_undirected;
532    goto &_edges_to if !wantarray;
533    $_[0]->[ _V ]->get_paths_by_ids([ &_edges_to ], &is_hyperedged && &is_directed);
534}
535
536sub successors {
537    my ($g, @args) = @_;
538    my ($V, $E) = @$g[ _V, _E ];
539    return if (my @i = $V->get_ids_by_paths(\@args)) != @args;
540    my @v = $E->successors(@i);
541    return @v if !wantarray;
542    map @$_, $V->get_paths_by_ids([ \@v ]);
543}
544
545sub predecessors {
546    goto &successors if &is_undirected;
547    my ($g, @args) = @_;
548    my ($V, $E) = @$g[ _V, _E ];
549    return if (my @i = $V->get_ids_by_paths(\@args)) != @args;
550    my @v = $E->predecessors(@i);
551    return @v if !wantarray;
552    map @$_, $V->get_paths_by_ids([ \@v ]);
553}
554
555sub _cessors_by_radius {
556    my ($radius, $method, $self_only_if_loop) = splice @_, -3, 3;
557    my ($g, @v) = @_;
558    require Set::Object;
559    my ($init, $next) = map Set::Object->new(@v), 1..2;
560    my $self = Set::Object->new(grep $g->has_edge($_, $_), @v) if $self_only_if_loop;
561    my ($got, $found) = map Set::Object->new, 1..2;
562    while (!defined $radius or $radius-- > 0) {
563	$found->insert($g->$method($next->members));
564	$next = $found->difference($got);
565	last if $next->is_null;  # Leave if no new found.
566	$got->insert($next->members);
567	$found->clear;
568    }
569    $got->remove($init->difference($self)->members) if $self_only_if_loop;
570    $got->${ wantarray ? \'members' : \'size' };
571}
572
573sub all_successors {
574    &expect_directed;
575    push @_, undef, 'successors', 0;
576    goto &_cessors_by_radius;
577}
578
579sub successors_by_radius {
580    &expect_directed;
581    push @_, 'successors', 0;
582    goto &_cessors_by_radius;
583}
584
585sub all_predecessors {
586    &expect_directed;
587    push @_, undef, 'predecessors', 0;
588    goto &_cessors_by_radius;
589}
590
591sub predecessors_by_radius {
592    &expect_directed;
593    push @_, 'predecessors', 0;
594    goto &_cessors_by_radius;
595}
596
597sub neighbours_by_radius {
598    push @_, 'neighbours', 1;
599    goto &_cessors_by_radius;
600}
601*neighbors_by_radius = \&neighbours_by_radius;
602
603sub neighbours {
604    require Set::Object;
605    my $s = Set::Object->new(&successors);
606    $s->insert(&predecessors) if &is_directed;
607    $s->${ wantarray ? \'members' : \'size' };
608}
609*neighbors = \&neighbours;
610
611sub all_neighbours {
612    push @_, undef, 'neighbours', 1;
613    goto &_cessors_by_radius;
614}
615*all_neighbors = \&all_neighbours;
616
617sub all_reachable {
618    &directed ? goto &all_successors : goto &all_neighbors;
619}
620
621sub reachable_by_radius {
622    &directed ? goto &successors_by_radius : goto &neighbors_by_radius;
623}
624
625sub delete_edge {
626    &expect_non_unionfind;
627    my $g = $_[0];
628    return $g if (my @i = &_vertex_ids) != @_ - 1;
629    @i = sort @i if &is_undirected;
630    return $g unless @i and $g->[ _E ]->del_path( \@i );
631    $g->[ _G ]++;
632    return $g;
633}
634
635sub delete_vertex {
636    &expect_non_unionfind;
637    my $g = $_[0];
638    return $g if @_ != 2;
639    my $V = $g->[ _V ];
640    return $g unless defined $V->has_path($_[1]);
641    # TODO: _edges_at is excruciatingly slow (rt.cpan.org 92427)
642    my $E = $g->[ _E ];
643    $E->del_path( $_ ) for &_edges_at;
644    $V->del_path($_[1]);
645    $g->[ _G ]++;
646    return $g;
647}
648
649sub get_vertex_count {
650    my $g = shift;
651    $g->[ _V ]->_get_path_count( @_ );
652}
653
654sub get_edge_count {
655    my $g = $_[0];
656    return 0 if (my @i = &_vertex_ids) != @_ - 1;
657    @i = sort @i if &is_undirected;
658    $g->[ _E ]->_get_path_count( \@i );
659}
660
661sub delete_vertices {
662    &expect_non_unionfind;
663    my $g = shift;
664    while (@_) {
665	my $v = shift @_;
666	$g->delete_vertex($v);
667    }
668    return $g;
669}
670
671sub delete_edges {
672    &expect_non_unionfind;
673    my $g = shift;
674    while (@_) {
675	my ($u, $v) = splice @_, 0, 2;
676	$g->delete_edge($u, $v);
677    }
678    return $g;
679}
680
681###
682# Degrees.
683#
684
685sub in_degree {
686    my $g = $_[0];
687    return undef unless @_ > 1 && &has_vertex;
688    my $in = 0;
689    $in += $g->get_edge_count( @$_ ) for &edges_to;
690    $in++ if &is_undirected and &is_self_loop_vertex;
691    return $in;
692}
693
694sub out_degree {
695    my $g = $_[0];
696    return undef unless @_ > 1 && &has_vertex;
697    my $out = 0;
698    $out += $g->get_edge_count( @$_ ) for &edges_from;
699    $out++ if &is_undirected and &is_self_loop_vertex;
700    return $out;
701}
702
703sub _total_degree {
704    return undef unless @_ > 1 && &has_vertex;
705    &is_undirected ? &in_degree : &in_degree - &out_degree;
706}
707
708sub degree {
709    goto &_total_degree if @_ > 1;
710    return 0 if &is_directed;
711    my $g = $_[0];
712    my $total = 0;
713    $total += $g->_total_degree( $_ ) for &_vertices05;
714    return $total;
715}
716
717*vertex_degree = \&degree;
718
719sub is_sink_vertex {
720    return 0 unless @_ > 1;
721    &successors == 0 && &predecessors > 0;
722}
723
724sub is_source_vertex {
725    return 0 unless @_ > 1;
726    &predecessors == 0 && &successors > 0;
727}
728
729sub is_successorless_vertex {
730    return 0 unless @_ > 1;
731    &successors == 0;
732}
733
734sub is_predecessorless_vertex {
735    return 0 unless @_ > 1;
736    &predecessors == 0;
737}
738
739sub is_successorful_vertex {
740    return 0 unless @_ > 1;
741    &successors > 0;
742}
743
744sub is_predecessorful_vertex {
745    return 0 unless @_ > 1;
746    &predecessors > 0;
747}
748
749sub is_isolated_vertex {
750    return 0 unless @_ > 1;
751    &predecessors == 0 && &successors == 0;
752}
753
754sub is_interior_vertex {
755    return 0 unless @_ > 1;
756    my $s = &successors;
757    $s-- if my $isl = &is_self_loop_vertex;
758    return 0 if $s == 0;
759    return $s > 0 if &is_undirected;
760    my $p = &predecessors;
761    $p-- if $isl;
762    $p > 0;
763}
764
765sub is_exterior_vertex {
766    return 0 unless @_ > 1;
767    &predecessors == 0 || &successors == 0;
768}
769
770sub is_self_loop_vertex {
771    return 0 unless @_ > 1;
772    return 1 if grep $_ eq $_[1], &successors; # @todo: multiedges
773    return 0;
774}
775
776for my $p (qw(
777    is_sink_vertex
778    is_source_vertex
779    is_successorless_vertex
780    is_predecessorless_vertex
781    is_successorful_vertex
782    is_predecessorful_vertex
783    is_isolated_vertex
784    is_interior_vertex
785    is_exterior_vertex
786    is_self_loop_vertex
787)) {
788    no strict 'refs';
789    (my $m = $p) =~ s/^is_(.*)ex$/${1}ices/;
790    *$m = sub { my $g = $_[0]; grep $g->$p($_), &_vertices05 };
791}
792
793###
794# Paths and cycles.
795#
796
797sub add_path {
798    my $g = shift;
799    my $u = shift;
800    my @edges;
801    while (@_) {
802	my $v = shift;
803	push @edges, [ $u, $v ];
804	$u = $v;
805    }
806    $g->add_edges(@edges);
807    return $g;
808}
809
810sub delete_path {
811    &expect_non_unionfind;
812    my $g = shift;
813    my $u = shift;
814    while (@_) {
815	my $v = shift;
816	$g->delete_edge($u, $v);
817	$u = $v;
818    }
819    return $g;
820}
821
822sub has_path {
823    my $g = shift;
824    my $u = shift;
825    while (@_) {
826	my $v = shift;
827	return 0 unless $g->has_edge($u, $v);
828	$u = $v;
829    }
830    return $g;
831}
832
833sub add_cycle {
834    push @_, $_[1];
835    goto &add_path;
836}
837
838sub delete_cycle {
839    &expect_non_unionfind;
840    push @_, $_[1];
841    goto &delete_path;
842}
843
844sub has_cycle {
845    return 0 if @_ == 1;
846    push @_, $_[1];
847    goto &has_path;
848}
849
850*has_this_cycle = \&has_cycle;
851
852sub has_a_cycle {
853    my $g = shift;
854    require Graph::Traversal::DFS;
855    my $t = Graph::Traversal::DFS->new($g, has_a_cycle => 1, @_);
856    $t->dfs;
857    return $t->get_state('has_a_cycle');
858}
859
860sub find_a_cycle {
861    require Graph::Traversal::DFS;
862    my @r = ( back_edge => \&Graph::Traversal::find_a_cycle);
863    push @r,
864      down_edge => \&Graph::Traversal::find_a_cycle
865	if &is_undirected;
866    my $g = shift;
867    my $t = Graph::Traversal::DFS->new($g, @r, @_);
868    $t->dfs;
869    $t->has_state('a_cycle') ? @{ $t->get_state('a_cycle') } : ();
870}
871
872###
873# Attributes.
874
875my @generic_methods = (
876    [ 'set_attribute', \&_set_attribute ],
877    [ 'set_attributes', \&_set_attributes ],
878    [ 'has_attributes', \&_has_attributes ],
879    [ 'has_attribute', \&_has_attribute ],
880    [ 'get_attributes', \&_get_attributes ],
881    [ 'get_attribute', \&_get_attribute ],
882    [ 'get_attribute_names', \&_get_attribute_names ],
883    [ 'get_attribute_values', \&_get_attribute_values ],
884    [ 'delete_attributes', \&_delete_attributes ],
885    [ 'delete_attribute', \&_delete_attribute ],
886);
887my %entity2offset = (vertex => _V, edge => _E);
888my %entity2args = (edge => '_vertex_ids');
889for my $entity (qw(vertex edge)) {
890    no strict 'refs';
891    my $expect_non = \&{ "expect_non_multi${entity}" };
892    my $expect_yes = \&{ "expect_multi${entity}" };
893    my $args_non = \&{ $entity2args{$entity} } if $entity2args{$entity};
894    my $args_yes = \&{ $entity2args{$entity}.'_multi' } if $entity2args{$entity};
895    my $offset = $entity2offset{$entity};
896    for my $t (@generic_methods) {
897	my ($raw, $func) = @$t;
898	my ($first, $rest) = ($raw =~ /^(\w+?)_(.+)/);
899	my $m = join '_', $first, $entity, $rest;
900	my $is_vertex = $entity eq 'vertex';
901	*$m = sub {
902	    &$expect_non; push @_, 0, $entity, $offset, $args_non, $is_vertex; goto &$func;
903	};
904	*{$m.'_by_id'} = sub {
905	    &$expect_yes; push @_, 1, $entity, $offset, $args_yes, $is_vertex; goto &$func;
906	};
907    }
908}
909
910sub _munge_args {
911    my ($is_vertex, $is_multi, $is_undirected, @args) = @_;
912    return \@args if !$is_vertex and !$is_undirected and !$is_multi;
913    return [ sort @args ] if !$is_vertex and !$is_multi;
914    return @args if $is_vertex;
915    my $id = pop @args;
916    ($is_undirected ? [ sort @args ] : \@args, $id);
917}
918
919sub _set_attribute {
920    my ($is_multi, $entity, $offset, $args, $is_vertex) = splice @_, -5, 5;
921    my $value = pop;
922    my $attr = pop;
923    no strict 'refs';
924    &{ 'add_' . $entity . ($is_multi ? '_by_id' : '') } unless &{ 'has_' . $entity . ($is_multi ? '_by_id' : '') };
925    my @args = ($entity eq 'edge') ? &$args : @_[1..$#_];
926    @args = _munge_args($is_vertex, $is_multi, &is_undirected, @args);
927    $_[0]->[ $offset ]->_set_path_attr( @args, $attr, $value );
928}
929
930sub _set_attributes {
931    my ($is_multi, $entity, $offset, $args, $is_vertex) = splice @_, -5, 5;
932    my $attr = pop;
933    no strict 'refs';
934    &{ 'add_' . $entity . ($is_multi ? '_by_id' : '') } unless &{ 'has_' . $entity . ($is_multi ? '_by_id' : '') };
935    my @args = ($entity eq 'edge') ? &$args : @_[1..$#_];
936    @args = _munge_args($is_vertex, $is_multi, &is_undirected, @args);
937    $_[0]->[ $offset ]->_set_path_attrs( @args, $attr );
938}
939
940sub _has_attributes {
941    my ($is_multi, $entity, $offset, $args, $is_vertex) = splice @_, -5, 5;
942    no strict 'refs';
943    return 0 unless &{ 'has_' . $entity . ($is_multi ? '_by_id' : '') };
944    my @args = ($entity eq 'edge') ? &$args : @_[1..$#_];
945    @args = _munge_args($is_vertex, $is_multi, &is_undirected, @args);
946    $_[0]->[ $offset ]->_has_path_attrs( @args );
947}
948
949sub _has_attribute {
950    my ($is_multi, $entity, $offset, $args, $is_vertex) = splice @_, -5, 5;
951    my $attr = pop;
952    no strict 'refs';
953    return 0 unless &{ 'has_' . $entity . ($is_multi ? '_by_id' : '') };
954    my @args = ($entity eq 'edge') ? &$args : @_[1..$#_];
955    @args = _munge_args($is_vertex, $is_multi, &is_undirected, @args);
956    $_[0]->[ $offset ]->_has_path_attr( @args, $attr );
957}
958
959sub _get_attributes {
960    my ($is_multi, $entity, $offset, $args, $is_vertex) = splice @_, -5, 5;
961    no strict 'refs';
962    return undef unless &{ 'has_' . $entity . ($is_multi ? '_by_id' : '') };
963    my @args = ($entity eq 'edge') ? &$args : @_[1..$#_];
964    @args = _munge_args($is_vertex, $is_multi, &is_undirected, @args);
965    scalar $_[0]->[ $offset ]->_get_path_attrs( @args );
966}
967
968sub _get_attribute {
969    my ($is_multi, $entity, $offset, $args, $is_vertex) = splice @_, -5, 5;
970    no strict 'refs';
971    my $attr = pop;
972    return undef unless &{ 'has_' . $entity . ($is_multi ? '_by_id' : '') };
973    my @args = ($entity eq 'edge') ? &$args : @_[1..$#_];
974    @args = _munge_args($is_vertex, $is_multi, &is_undirected, @args);
975    scalar $_[0]->[ $offset ]->_get_path_attr( @args, $attr );
976}
977
978sub _get_attribute_names {
979    my ($is_multi, $entity, $offset, $args, $is_vertex) = splice @_, -5, 5;
980    no strict 'refs';
981    return unless &{ 'has_' . $entity . ($is_multi ? '_by_id' : '') };
982    my @args = ($entity eq 'edge') ? &$args : @_[1..$#_];
983    @args = _munge_args($is_vertex, $is_multi, &is_undirected, @args);
984    $_[0]->[ $offset ]->_get_path_attr_names( @args );
985}
986
987sub _get_attribute_values {
988    my ($is_multi, $entity, $offset, $args, $is_vertex) = splice @_, -5, 5;
989    no strict 'refs';
990    return unless &{ 'has_' . $entity . ($is_multi ? '_by_id' : '') };
991    my @args = ($entity eq 'edge') ? &$args : @_[1..$#_];
992    @args = _munge_args($is_vertex, $is_multi, &is_undirected, @args);
993    $_[0]->[ $offset ]->_get_path_attr_values( @args );
994}
995
996sub _delete_attributes {
997    my ($is_multi, $entity, $offset, $args, $is_vertex) = splice @_, -5, 5;
998    no strict 'refs';
999    return undef unless &{ 'has_' . $entity . ($is_multi ? '_by_id' : '') };
1000    my @args = ($entity eq 'edge') ? &$args : @_[1..$#_];
1001    @args = _munge_args($is_vertex, $is_multi, &is_undirected, @args);
1002    $_[0]->[ $offset ]->_del_path_attrs( @args );
1003}
1004
1005sub _delete_attribute {
1006    my ($is_multi, $entity, $offset, $args, $is_vertex) = splice @_, -5, 5;
1007    my $attr = pop;
1008    no strict 'refs';
1009    return undef unless &{ 'has_' . $entity . ($is_multi ? '_by_id' : '') };
1010    my @args = ($entity eq 'edge') ? &$args : @_[1..$#_];
1011    @args = _munge_args($is_vertex, $is_multi, &is_undirected, @args);
1012    $_[0]->[ $offset ]->_del_path_attr( @args, $attr );
1013}
1014
1015sub add_vertices {
1016    my ($g, @v) = @_;
1017    if (&is_multivertexed) {
1018	$g->add_vertex_by_id($_, _GEN_ID) for @v;
1019	return $g;
1020    }
1021    my @i = $g->[ _V ]->set_paths(@v);
1022    $g->[ _G ]++;
1023    return $g if !&has_union_find;
1024    $g->[ _U ]->add(@i);
1025    $g;
1026}
1027
1028sub add_edges {
1029    my ($g, @args) = @_;
1030    my @edges;
1031    while (defined(my $u = shift @args)) {
1032	push @edges, ref $u eq 'ARRAY' ? $u : @args ? [ $u, shift @args ]
1033	    : __carp_confess "Graph::add_edges: missing end vertex";
1034    }
1035    if (&is_multiedged) {
1036	$g->add_edge_by_id(@$_, _GEN_ID) for @edges;
1037	return $g;
1038    }
1039    my $uf = &has_union_find;
1040    my $deep = &is_hyperedged && &is_directed;
1041    my @paths = $g->[ _V ]->get_ids_by_paths(\@edges, 1, 1 + ($deep ? 1 : 0));
1042    @paths = map [ sort @$_ ], @paths if &is_undirected;
1043    $g->[ _E ]->set_paths( @paths );
1044    $uf->union(@paths) if $uf;
1045    $g->[ _G ]++;
1046    return $g;
1047}
1048
1049sub rename_vertex {
1050    my $g = shift;
1051    $g->[ _V ]->rename_path(@_);
1052    return $g;
1053}
1054
1055sub rename_vertices {
1056    my ($g, $code) = @_;
1057    my %seen;
1058    $g->rename_vertex($_, $code->($_))
1059	for grep !$seen{$_}++, $g->[ _V ]->paths;
1060    return $g;
1061}
1062
1063sub as_hashes {
1064    my ($g) = @_;
1065    my (%v, %e, @e);
1066    my ($is_hyper, $is_directed)= (&is_hyperedged, &is_directed);
1067    if (&is_multivertexed) {
1068        for my $v ($g->unique_vertices) {
1069            $v{$v} = {
1070                map +($_ => $g->get_vertex_attributes_by_id($v, $_) || {}),
1071                    $g->get_multivertex_ids($v)
1072            };
1073        }
1074    } else {
1075        %v = map +($_ => $g->get_vertex_attributes($_) || {}), $g->unique_vertices;
1076    }
1077    my $multi_e = &is_multiedged;
1078    for my $e ($g->edges) {
1079	my $edge_attr = {
1080	    $multi_e
1081		? map +($_ => $g->get_edge_attributes_by_id(@$e, $_) || {}),
1082		    $g->get_multiedge_ids(@$e)
1083		: %{ $g->get_edge_attributes(@$e)||{} }
1084	};
1085	if ($is_hyper) {
1086	    my %h = (attributes => $edge_attr);
1087	    if ($is_directed) {
1088		@h{qw(predecessors successors)} = @$e;
1089	    } else {
1090		$h{vertices} = $e;
1091	    }
1092	    push @e, \%h;
1093	} else {
1094	    $e{ $e->[0] }{ $e->[1] } = $edge_attr;
1095	    $e{ $e->[1] }{ $e->[0] } = $edge_attr if !$is_directed;
1096	}
1097    }
1098    ( \%v, $is_hyper ? \@e : \%e );
1099}
1100
1101sub ingest {
1102    my ($g, $g2) = @_;
1103    for my $v ($g2->vertices) {
1104        if (&is_multivertexed) {
1105            $g->set_vertex_attributes_by_id($v, $_, $g2->get_vertex_attributes_by_id($v, $_))
1106                for $g2->get_multivertex_ids($v);
1107        } else {
1108            $g->set_vertex_attributes($v, $g2->get_vertex_attributes($v));
1109        }
1110        if (&is_multiedged) {
1111            for my $e ($g2->edges_from($v)) {
1112                $g->set_edge_attributes_by_id(@$e, $_, $g2->get_edge_attributes_by_id(@$e, $_))
1113                    for $g2->get_multiedge_ids(@$e);
1114            }
1115        } else {
1116            $g->set_edge_attributes(@$_, $g2->get_edge_attributes(@$_))
1117                for $g2->edges_from($v);
1118        }
1119    }
1120    $g;
1121}
1122
1123###
1124# More constructors.
1125#
1126
1127sub copy {
1128    my ($g, @args) = @_;
1129    my %opt = _get_options( \@args );
1130    no strict 'refs';
1131    my $c = (ref $g)->new(map +($_ => &$_ ? 1 : 0), @GRAPH_PROPS_COPIED);
1132    $c->add_vertices(&isolated_vertices);
1133    $c->add_edges(&_edges05);
1134    return $c;
1135}
1136
1137*copy_graph = \&copy;
1138
1139sub _deep_copy_best {
1140    _can_deep_copy_Storable()
1141        ? _deep_copy_Storable(@_) : _deep_copy_DataDumper(@_);
1142}
1143
1144sub _deep_copy_Storable {
1145    my $g = shift;
1146    require Safe;   # For deep_copy().
1147    my $safe = Safe->new;
1148    $safe->permit(qw/:load/);
1149    local $Storable::Deparse = 1;
1150    local $Storable::Eval = sub { $safe->reval($_[0]) };
1151    return Storable::thaw(Storable::freeze($g));
1152}
1153
1154sub _deep_copy_DataDumper {
1155    my $g = shift;
1156    require Data::Dumper;
1157    my $d = Data::Dumper->new([$g]);
1158    use vars qw($VAR1);
1159    $d->Purity(1)->Terse(1)->Deepcopy(1);
1160    $d->Deparse(1) if $] >= 5.008;
1161    eval $d->Dump;
1162}
1163
1164sub deep_copy {
1165    local $. = $.;
1166    my $g2 = _deep_copy_best(@_);
1167    $g2->[ _V ]->reindex if grep ref, &_vertices05;
1168    $g2;
1169}
1170
1171*deep_copy_graph = \&deep_copy;
1172
1173sub transpose_edge {
1174    my $g = $_[0];
1175    return $g if !&is_directed;
1176    return undef unless &has_edge;
1177    my $c = &get_edge_count;
1178    my $a = &get_edge_attributes;
1179    my @e = reverse @_[1..$#_];
1180    &delete_edge unless $g->has_edge( @e );
1181    $g->add_edges(map \@e, 1..$c);
1182    $g->set_edge_attributes(@e, $a) if $a;
1183    return $g;
1184}
1185
1186sub transpose_graph {
1187    my $t = &copy;
1188    return $t if !&directed;
1189    $t->transpose_edge(@$_) for &_edges05;
1190    return $t;
1191}
1192
1193*transpose = \&transpose_graph;
1194
1195sub complete_graph {
1196    my $directed = &is_directed;
1197    my $c = &new;
1198    my @v = &_vertices05;
1199    my @edges;
1200    for (my $i = $#v; $i >= 0; $i-- ) {
1201	push @edges, map +([$v[$i], $v[$_]], $directed ? [$v[$_], $v[$i]] : ()),
1202	    0..$i - 1;
1203    }
1204    $c->add_edges(@edges);
1205    return $c;
1206}
1207
1208*complement = \&complement_graph;
1209
1210sub complement_graph {
1211    my $c = &complete_graph;
1212    $c->delete_edge(@$_) for &edges;
1213    return $c;
1214}
1215
1216*complete = \&complete_graph;
1217
1218sub subgraph {
1219  my ($g, $src, $dst) = @_;
1220  __carp_confess "Graph::subgraph: need src and dst array references"
1221    unless ref $src eq 'ARRAY' && (!defined($dst) or ref $dst eq 'ARRAY');
1222  require Set::Object;
1223  my $s = $g->new;
1224  my @u = grep $g->has_vertex($_), @$src;
1225  my $v = Set::Object->new($dst ? grep $g->has_vertex($_), @$dst : @u);
1226  $s->add_vertices(@u, $dst ? $v->members : ());
1227  my $directed = &is_directed;
1228  $s->add_edges(grep $v->contains($directed ? $_->[1] : @$_), $g->edges_from(@u));
1229  return $s;
1230}
1231
1232###
1233# Transitivity.
1234#
1235
1236sub is_transitive {
1237    my $g = shift;
1238    require Graph::TransitiveClosure;
1239    Graph::TransitiveClosure::is_transitive($g);
1240}
1241
1242###
1243# Weighted vertices.
1244#
1245
1246my $defattr = 'weight';
1247
1248sub _defattr {
1249    return $defattr;
1250}
1251
1252sub add_weighted_vertex {
1253    &expect_non_multivertexed;
1254    push @_, $defattr, pop;
1255    goto &set_vertex_attribute;
1256}
1257
1258sub add_weighted_vertices {
1259    &expect_non_multivertexed;
1260    my $g = shift;
1261    while (@_) {
1262	my ($v, $w) = splice @_, 0, 2;
1263	$g->set_vertex_attribute($v, $defattr, $w);
1264    }
1265}
1266
1267sub get_vertex_weight {
1268    &expect_non_multivertexed;
1269    push @_, $defattr;
1270    goto &get_vertex_attribute;
1271}
1272
1273sub has_vertex_weight {
1274    &expect_non_multivertexed;
1275    push @_, $defattr;
1276    goto &has_vertex_attribute;
1277}
1278
1279sub set_vertex_weight {
1280    &expect_non_multivertexed;
1281    push @_, $defattr, pop;
1282    goto &set_vertex_attribute;
1283}
1284
1285sub delete_vertex_weight {
1286    &expect_non_multivertexed;
1287    push @_, $defattr;
1288    goto &delete_vertex_attribute;
1289}
1290
1291sub add_weighted_vertex_by_id {
1292    &expect_multivertexed;
1293    push @_, $defattr, pop;
1294    goto &set_vertex_attribute_by_id;
1295}
1296
1297sub add_weighted_vertices_by_id {
1298    &expect_multivertexed;
1299    my $g = shift;
1300    my $id = pop;
1301    while (@_) {
1302	my ($v, $w) = splice @_, 0, 2;
1303	$g->add_vertex_by_id($v, $id);
1304	$g->set_vertex_attribute_by_id($v, $id, $defattr, $w);
1305    }
1306}
1307
1308sub get_vertex_weight_by_id {
1309    &expect_multivertexed;
1310    push @_, $defattr;
1311    goto &get_vertex_attribute_by_id;
1312}
1313
1314sub has_vertex_weight_by_id {
1315    &expect_multivertexed;
1316    push @_, $defattr;
1317    goto &has_vertex_attribute_by_id;
1318}
1319
1320sub set_vertex_weight_by_id {
1321    &expect_multivertexed;
1322    push @_, $defattr, pop;
1323    goto &set_vertex_attribute_by_id;
1324}
1325
1326sub delete_vertex_weight_by_id {
1327    &expect_multivertexed;
1328    push @_, $defattr;
1329    goto &delete_vertex_attribute_by_id;
1330}
1331
1332###
1333# Weighted edges.
1334#
1335
1336sub add_weighted_edge {
1337    &expect_non_multiedged;
1338    push @_, $defattr, pop;
1339    goto &set_edge_attribute;
1340}
1341
1342sub add_weighted_edges {
1343    &expect_non_multiedged;
1344    my $g = shift;
1345    while (@_) {
1346	my ($u, $v, $w) = splice @_, 0, 3;
1347	$g->set_edge_attribute($u, $v, $defattr, $w);
1348    }
1349}
1350
1351sub add_weighted_edges_by_id {
1352    &expect_multiedged;
1353    my $g = shift;
1354    my $id = pop;
1355    while (@_) {
1356	my ($u, $v, $w) = splice @_, 0, 3;
1357	$g->set_edge_attribute_by_id($u, $v, $id, $defattr, $w);
1358    }
1359}
1360
1361sub add_weighted_path {
1362    &expect_non_multiedged;
1363    my $g = shift;
1364    my $u = shift;
1365    while (@_) {
1366	my ($w, $v) = splice @_, 0, 2;
1367	$g->set_edge_attribute($u, $v, $defattr, $w);
1368	$u = $v;
1369    }
1370}
1371
1372sub get_edge_weight {
1373    &expect_non_multiedged;
1374    push @_, $defattr;
1375    goto &get_edge_attribute;
1376}
1377
1378sub has_edge_weight {
1379    &expect_non_multiedged;
1380    push @_, $defattr;
1381    goto &has_edge_attribute;
1382}
1383
1384sub set_edge_weight {
1385    &expect_non_multiedged;
1386    push @_, $defattr, pop;
1387    goto &set_edge_attribute;
1388}
1389
1390sub delete_edge_weight {
1391    &expect_non_multiedged;
1392    push @_, $defattr;
1393    goto &delete_edge_attribute;
1394}
1395
1396sub add_weighted_edge_by_id {
1397    &expect_multiedged;
1398    push @_, $defattr, pop;
1399    goto &set_edge_attribute_by_id;
1400}
1401
1402sub add_weighted_path_by_id {
1403    &expect_multiedged;
1404    my $g = shift;
1405    my $id = pop;
1406    my $u = shift;
1407    while (@_) {
1408	my ($w, $v) = splice @_, 0, 2;
1409	$g->set_edge_attribute_by_id($u, $v, $id, $defattr, $w);
1410	$u = $v;
1411    }
1412}
1413
1414sub get_edge_weight_by_id {
1415    &expect_multiedged;
1416    push @_, $defattr;
1417    goto &get_edge_attribute_by_id;
1418}
1419
1420sub has_edge_weight_by_id {
1421    &expect_multiedged;
1422    push @_, $defattr;
1423    goto &has_edge_attribute_by_id;
1424}
1425
1426sub set_edge_weight_by_id {
1427    &expect_multiedged;
1428    push @_, $defattr, pop;
1429    goto &set_edge_attribute_by_id;
1430}
1431
1432sub delete_edge_weight_by_id {
1433    &expect_multiedged;
1434    push @_, $defattr;
1435    goto &delete_edge_attribute_by_id;
1436}
1437
1438###
1439# Error helpers.
1440#
1441
1442my %expected;
1443@expected{qw(directed undirected acyclic)} = qw(undirected directed cyclic);
1444
1445sub _expected {
1446    my $exp = shift;
1447    my $got = @_ ? shift : $expected{$exp};
1448    $got = defined $got ? ", got $got" : "";
1449    if (my @caller2 = caller(2)) {
1450	die "$caller2[3]: expected $exp graph$got, at $caller2[1] line $caller2[2].\n";
1451    } else {
1452	my @caller1 = caller(1); # uncoverable statement
1453	die "$caller1[3]: expected $exp graph$got, at $caller1[1] line $caller1[2].\n"; # uncoverable statement
1454    }
1455}
1456
1457sub expect_no_args {
1458    my $g = shift;
1459    return unless @_;
1460    my @caller1 = caller(1); # uncoverable statement
1461    die "$caller1[3]: expected no arguments, got " . scalar @_ . ", at $caller1[1] line $caller1[2]\n"; # uncoverable statement
1462}
1463
1464sub expect_undirected {
1465    _expected('undirected') unless &is_undirected;
1466}
1467
1468sub expect_directed {
1469    _expected('directed') unless &is_directed;
1470}
1471
1472sub expect_acyclic {
1473    _expected('acyclic') unless &is_acyclic;
1474}
1475
1476sub expect_dag {
1477    my @got;
1478    push @got, 'undirected' unless &is_directed;
1479    push @got, 'cyclic'     unless &is_acyclic;
1480    _expected('directed acyclic', "@got") if @got;
1481}
1482
1483sub expect_hyperedged {
1484    _expected('hyperedged') unless &is_hyperedged;
1485}
1486
1487sub expect_multivertexed {
1488    _expected('multivertexed') unless &is_multivertexed;
1489}
1490*expect_multivertex = \&expect_multivertexed;
1491
1492sub expect_non_multivertexed {
1493    _expected('non-multivertexed') if &is_multivertexed;
1494}
1495*expect_non_multivertex = \&expect_non_multivertexed;
1496
1497sub expect_non_multiedged {
1498    _expected('non-multiedged') if &is_multiedged;
1499}
1500*expect_non_multiedge = \&expect_non_multiedged;
1501
1502sub expect_multiedged {
1503    _expected('multiedged') unless &is_multiedged;
1504}
1505*expect_multiedge = \&expect_multiedged;
1506
1507sub expect_non_unionfind {
1508    _expected('non-unionfind') if &has_union_find;
1509}
1510
1511sub _get_options {
1512    my @caller = caller(1);
1513    unless (@_ == 1 && ref $_[0] eq 'ARRAY') {
1514	die "$caller[3]: internal error: should be called with only one array ref argument, at $caller[1] line $caller[2].\n";
1515    }
1516    my @opt = @{ $_[0] };
1517    unless (@opt  % 2 == 0) {
1518	die "$caller[3]: expected an options hash, got a non-even number of arguments, at $caller[1] line $caller[2].\n"; # uncoverable statement
1519    }
1520    return @opt;
1521}
1522
1523###
1524# Random constructors and accessors.
1525#
1526
1527sub __fisher_yates_shuffle (@) {
1528    # From perlfaq4, but modified to be non-modifying.
1529    my @a = @_;
1530    my $i = @a;
1531    while ($i--) {
1532	my $j = int rand ($i+1);
1533	@a[$i,$j] = @a[$j,$i];
1534    }
1535    return @a;
1536}
1537
1538BEGIN {
1539    sub _shuffle(@);
1540    # Workaround for the Perl bug [perl #32383] where -d:Dprof and
1541    # List::Util::shuffle do not like each other: if any debugging
1542    # (-d) flags are on, fall back to our own Fisher-Yates shuffle.
1543    # The bug was fixed by perl changes #26054 and #26062, which
1544    # went to Perl 5.9.3.  If someone tests this with a pre-5.9.3
1545    # bleadperl that calls itself 5.9.3 but doesn't yet have the
1546    # patches, oh, well.
1547    *_shuffle = $^P && $] < 5.009003 ?
1548	\&__fisher_yates_shuffle : do { require List::Util; \&List::Util::shuffle };
1549}
1550
1551sub random_graph {
1552    my $class = (@_ % 2) == 0 ? 'Graph' : shift;
1553    my %opt = _get_options( \@_ );
1554    __carp_confess "Graph::random_graph: argument 'vertices' missing or undef"
1555	unless defined $opt{vertices};
1556    srand delete $opt{random_seed} if exists $opt{random_seed};
1557    my $random_edge = delete $opt{random_edge} if exists $opt{random_edge};
1558    my @V;
1559    if (my $ref = ref $opt{vertices}) {
1560	__carp_confess "Graph::random_graph: argument 'vertices' illegal"
1561	    if $ref ne 'ARRAY';
1562	@V = @{ $opt{vertices} };
1563    } else {
1564	@V = 0..($opt{vertices} - 1);
1565    }
1566    delete $opt{vertices};
1567    my $V = @V;
1568    my $C = $V * ($V - 1) / 2;
1569    my $E;
1570    __carp_confess "Graph::random_graph: both arguments 'edges' and 'edges_fill' specified"
1571	if exists $opt{edges} && exists $opt{edges_fill};
1572    $E = exists $opt{edges_fill} ? $opt{edges_fill} * $C : $opt{edges};
1573    delete $opt{edges};
1574    delete $opt{edges_fill};
1575    my $g = $class->new(%opt);
1576    $g->add_vertices(@V);
1577    return $g if $V < 2;
1578    $C *= 2 if my $is_directed = $g->directed;
1579    $E = $C / 2 unless defined $E;
1580    $E = int($E + 0.5);
1581    my $p = $E / $C;
1582    $random_edge = sub { $p } unless defined $random_edge;
1583    # print "V = $V, E = $E, C = $C, p = $p\n";
1584    __carp_confess "Graph::random_graph: needs to be countedged or multiedged ($E > $C)"
1585	if $p > 1.0 && !($g->countedged || $g->multiedged);
1586    # Shuffle the vertex lists so that the pairs at
1587    # the beginning of the lists are not more likely.
1588    my (%v1_v2, @edges);
1589    my @V1 = _shuffle @V;
1590    my @V2 = _shuffle @V;
1591 LOOP:
1592    while ($E) {
1593	for my $v1 (@V1) {
1594	    for my $v2 (@V2) {
1595		next if $v1 eq $v2; # TODO: allow self-loops?
1596		my $q = $random_edge->($g, $v1, $v2, $p);
1597		if ($q && ($q == 1 || rand() <= $q) &&
1598		    !exists $v1_v2{$v1}{$v2} &&
1599		    ($is_directed ? 1 : !exists $v1_v2{$v2}{$v1})) {
1600		    $v1_v2{$v1}{$v2} = undef;
1601		    push @edges, [ $v1, $v2 ];
1602		    $E--;
1603		    last LOOP unless $E;
1604		}
1605	    }
1606	}
1607    }
1608    $g->add_edges(@edges);
1609}
1610
1611sub random_vertex {
1612    my @V = &_vertices05;
1613    @V[rand @V];
1614}
1615
1616sub random_edge {
1617    my @E = &_edges05;
1618    @E[rand @E];
1619}
1620
1621sub random_successor {
1622    my @S = &successors;
1623    @S[rand @S];
1624}
1625
1626sub random_predecessor {
1627    my @P = &predecessors;
1628    @P[rand @P];
1629}
1630
1631###
1632# Algorithms.
1633#
1634
1635my $MST_comparator = sub { ($_[0] || 0) <=> ($_[1] || 0) };
1636
1637sub _MST_attr {
1638    my $attr = shift;
1639    my $attribute =
1640	exists $attr->{attribute}  ?
1641	    $attr->{attribute}  : $defattr;
1642    my $comparator =
1643	exists $attr->{comparator} ?
1644	    $attr->{comparator} : $MST_comparator;
1645    return ($attribute, $comparator);
1646}
1647
1648sub _MST_edges {
1649    my ($g, $attr) = @_;
1650    my ($attribute, $comparator) = _MST_attr($attr);
1651    map $_->[1],
1652        sort { $comparator->($a->[0], $b->[0], $a->[1], $b->[1]) }
1653             map [ $g->get_edge_attribute(@$_, $attribute), $_ ],
1654                 &_edges05;
1655}
1656
1657sub MST_Kruskal {
1658    &expect_undirected;
1659    my ($g, %attr) = @_;
1660    require Graph::UnionFind;
1661
1662    my $MST = Graph->new(directed => 0);
1663
1664    my $UF  = Graph::UnionFind->new;
1665    $UF->add(&_vertices05);
1666
1667    my @edges;
1668    for my $e ($g->_MST_edges(\%attr)) {
1669	my ($u, $v) = @$e; # TODO: hyperedges
1670	next if $UF->same( @$e );
1671	$UF->union([$u, $v]);
1672	push @edges, [ $u, $v ];
1673    }
1674    $MST->add_edges(@edges);
1675
1676    return $MST;
1677}
1678
1679sub _MST_add {
1680    my ($g, $h, $HF, $r, $attr, $unseen) = @_;
1681    $HF->add( Graph::MSTHeapElem->new( $r, $_, $g->get_edge_attribute( $r, $_, $attr ) ) )
1682	for grep exists $unseen->{ $_ }, $g->successors( $r );
1683}
1684
1685sub _next_alphabetic { shift; (sort               keys %{ $_[0] })[0] }
1686sub _next_numeric    { shift; (sort { $a <=> $b } keys %{ $_[0] })[0] }
1687sub _next_random     { shift; (values %{ $_[0] })[ rand keys %{ $_[0] } ] }
1688
1689sub _root_opt {
1690    my ($g, @args) = @_;
1691    my %opt = @args == 1 ? ( first_root => $args[0] ) : _get_options( \@args );
1692    my %unseen;
1693    my @unseen = $g->_vertices05;
1694    @unseen{ @unseen } = @unseen;
1695    @unseen = _shuffle @unseen;
1696    my $r;
1697    if (exists $opt{ start }) {
1698	$opt{ first_root } = delete $opt{ start };
1699	$opt{ next_root  } = undef;
1700    }
1701    if (exists $opt{ first_root }) {
1702	if (ref $opt{ first_root } eq 'CODE') {
1703	    $r = $opt{ first_root }->( $g, \%unseen );
1704	} else {
1705	    $r = $opt{ first_root };
1706	}
1707    } else {
1708	$r = shift @unseen;
1709    }
1710    my $next =
1711	exists $opt{ next_root } ?
1712	    $opt{ next_root } :
1713              $opt{ next_alphabetic } ?
1714                \&_next_alphabetic :
1715                  $opt{ next_numeric } ?
1716                    \&_next_numeric :
1717                      \&_next_random;
1718    my $code = ref $next eq 'CODE';
1719    my $attr = exists $opt{ attribute } ? $opt{ attribute } : $defattr;
1720    return ( \%opt, \%unseen, \@unseen, $r, $next, $code, $attr );
1721}
1722
1723sub _heap_walk {
1724    my ($g, $h, $add, $etc, $opt, $unseenh, $unseena, $r, $next, $code, $attr) = @_;
1725    require Heap::Fibonacci;
1726    my $HF = Heap::Fibonacci->new;
1727    while (defined $r) {
1728        # print "r = $r\n";
1729	$add->($g, $h, $HF, $r, $attr, $unseenh, $etc);
1730	delete $unseenh->{ $r };
1731	while (defined $HF->top) {
1732	    my $t = $HF->extract_top;
1733	    # use Data::Dumper; print "t = ", Dumper($t);
1734	    if (defined $t) {
1735		my ($u, $v, $w) = $t->val;
1736		# print "extracted top: $u $v $w\n";
1737		if (exists $unseenh->{ $v }) {
1738		    $h->set_edge_attribute($u, $v, $attr, $w);
1739		    delete $unseenh->{ $v };
1740		    $add->($g, $h, $HF, $v, $attr, $unseenh, $etc);
1741		}
1742	    }
1743	}
1744	return $h unless defined $next;
1745	$r = $code ? $next->( $g, $unseenh ) : shift @$unseena;
1746        last unless defined $r;
1747    }
1748    return $h;
1749}
1750
1751sub MST_Prim {
1752    &expect_undirected;
1753    require Graph::MSTHeapElem;
1754    $_[0]->_heap_walk(Graph->new(directed => 0), \&_MST_add, undef, &_root_opt);
1755}
1756
1757*MST_Dijkstra = \&MST_Prim;
1758
1759*minimum_spanning_tree = \&MST_Prim;
1760
1761###
1762# Cycle detection.
1763#
1764
1765*is_cyclic = \&has_a_cycle;
1766
1767sub is_acyclic {
1768    !&is_cyclic;
1769}
1770
1771sub is_dag {
1772    &is_directed && &is_acyclic ? 1 : 0;
1773}
1774
1775*is_directed_acyclic_graph = \&is_dag;
1776
1777###
1778# Simple DFS uses.
1779#
1780
1781sub topological_sort {
1782    my $g = shift;
1783    my %opt = _get_options( \@_ );
1784    my $eic = delete $opt{ empty_if_cyclic };
1785    my $hac;
1786    if ($eic) {
1787	$hac = $g->has_a_cycle;
1788    } else {
1789	$g->expect_dag;
1790    }
1791    require Graph::Traversal::DFS;
1792    my $t = Graph::Traversal::DFS->new($g, %opt);
1793    my @s = $t->dfs;
1794    $hac ? () : reverse @s;
1795}
1796
1797*toposort = \&topological_sort;
1798
1799sub _undirected_copy_compute {
1800  Graph->new(directed => 0, vertices => [&isolated_vertices], edges => [&_edges05]);
1801}
1802
1803sub undirected_copy {
1804    &expect_directed;
1805    return _check_cache($_[0], 'undirected_copy', [], \&_undirected_copy_compute);
1806}
1807
1808*undirected_copy_graph = \&undirected_copy;
1809
1810sub directed_copy {
1811    &expect_undirected;
1812    my @edges = &_edges05;
1813    Graph->new(directed => 1, vertices => [&isolated_vertices],
1814	edges => [@edges, map [reverse @$_], @edges]);
1815}
1816
1817*directed_copy_graph = \&directed_copy;
1818
1819###
1820# Cache or not.
1821#
1822
1823my %_cache_type =
1824    (
1825     'connectivity'        => ['_ccc'],
1826     'strong_connectivity' => ['_scc'],
1827     'biconnectivity'      => ['_bcc'],
1828     'SPT_Dijkstra'        => ['_spt_di', 'SPT_Dijkstra_root'],
1829     'SPT_Bellman_Ford'    => ['_spt_bf', 'SPT_Bellman_Ford_root'],
1830     'undirected_copy'     => ['_undirected'],
1831     'transitive_closure_matrix' => ['_tcm'],
1832    );
1833
1834for my $t (keys %_cache_type) {
1835    no strict 'refs';
1836    my @attr = @{ $_cache_type{$t} };
1837    *{$t."_clear_cache"} = sub { $_[0]->delete_graph_attribute($_) for @attr };
1838}
1839
1840sub _check_cache {
1841    my ($g, $type, $extra_vals, $code, @args) = @_;
1842    my $c = $_cache_type{$type};
1843    __carp_confess "Graph: unknown cache type '$type'" if !defined $c;
1844    my ($main_key, @extra_keys) = @$c;
1845    __carp_confess "Graph: wrong number of extra values (@extra_keys) vs (@$extra_vals)" if @extra_keys != @$extra_vals;
1846    my $a = $g->get_graph_attribute($main_key);
1847    __carp_confess "$c attribute set to unexpected value $a"
1848	if defined $a and ref $a ne 'ARRAY';
1849    unless (defined $a && $a->[ 0 ] == $g->[ _G ]) {
1850	$g->set_graph_attribute($main_key, $a = [ $g->[ _G ], $code->( $g, @args ) ]);
1851    }
1852    my $i = -1;
1853    my $extra_invalid = grep {
1854	my $v = $a->[ 1 ]->get_graph_attribute($_);
1855	$i++; # here so still incremented even if short-cut
1856	!defined $v or $v ne $extra_vals->[$i];
1857    } @extra_keys;
1858    if ($extra_invalid) {
1859	$g->set_graph_attribute($main_key, $a = [ $g->[ _G ], $code->( $g, @args ) ]);
1860    }
1861    return $a->[ 1 ];
1862}
1863
1864###
1865# Connected components.
1866#
1867
1868sub _connected_components_compute {
1869    my $g = $_[0];
1870    my %v2c;
1871    my @c;
1872    return [ [], {} ] unless my @v = $g->unique_vertices;
1873    if (my $UF = &has_union_find) {
1874	my $V  = $g->[ _V ];
1875	my @ids = $V->get_ids_by_paths(\@v, 0);
1876	my ($counter, %cc2counter) = 0;
1877	my @cc = $UF->find(@ids);
1878	for (my $i = 0; $i <= $#v; $i++) {
1879	    my $cc = $cc[$i];
1880	    __carp_confess "connected_component union-find did not have vertex '$v[$i]', please report"
1881		if !defined $cc;
1882	    $cc2counter{$cc} = $counter++ if !exists $cc2counter{$cc};
1883	    my $ci = $cc2counter{$cc};
1884	    $v2c{ $v[$i] } = $ci;
1885	    push @{ $c[$ci] }, $v[$i];
1886	}
1887    } else {
1888	require Graph::Traversal::DFS;
1889	my %r; @r{ @v } = @v;
1890	@c = [];
1891	my $t = Graph::Traversal::DFS->new(
1892	    $g,
1893	    first_root => sub { (each %r)[1] },
1894	    next_root  => sub { push @c, [] if keys %r; (each %r)[1]; },
1895	    pre => sub {
1896		my ($v, $t) = @_;
1897		$v2c{ $v } = $#c;
1898		push @{ $c[-1] }, $v;
1899		delete $r{ $v };
1900	    },
1901	    @_[1..$#_]
1902	);
1903	$t->dfs;
1904    }
1905    return [ \@c, \%v2c ];
1906}
1907
1908sub _connected_components {
1909    my $ccc = _check_cache($_[0], 'connectivity', [],
1910			   \&_connected_components_compute);
1911    return @{ $ccc };
1912}
1913
1914sub connected_component_by_vertex {
1915    &expect_undirected;
1916    (&_connected_components)[1]->{ $_[1] };
1917}
1918
1919sub connected_component_by_index {
1920    &expect_undirected;
1921    my $value = (&_connected_components)[0]->[$_[1]];
1922    $value ? @{ $value || _empty_array } : ();
1923}
1924
1925sub connected_components {
1926    &expect_undirected;
1927    @{ (&_connected_components)[0] };
1928}
1929
1930sub same_connected_components {
1931    &expect_undirected;
1932    my ($g, @args) = @_;
1933    my @components;
1934    if (my $UF = &has_union_find) {
1935	my @ids = &_vertex_ids;
1936	return 0 if @ids != @args;
1937	@components = $UF->find(@ids);
1938    } else {
1939	@components = @{ (&_connected_components)[1] }{ @args };
1940    }
1941    return 0 if grep !defined, @components;
1942    require List::Util;
1943    List::Util::uniq( @components ) == 1;
1944}
1945
1946sub _super_component { join("+", sort @_) }
1947
1948sub connected_graph {
1949    &expect_undirected;
1950    my ($g, %opt) = @_;
1951    my $cg = Graph->new(undirected => 1);
1952    if ($g->has_union_find && $g->vertices == 1) {
1953	# TODO: super_component?
1954	$cg->add_vertices($g->vertices);
1955    } else {
1956	my $sc_cb = $opt{super_component} || \&_super_component;
1957	$cg->set_vertex_attribute(scalar $sc_cb->(@$_), 'subvertices', $_)
1958	    for $g->connected_components;
1959    }
1960    return $cg;
1961}
1962
1963sub is_connected {
1964    &expect_undirected;
1965    return @{ (&_connected_components)[0] } == 1;
1966}
1967
1968sub is_weakly_connected {
1969    &expect_directed;
1970    splice @_, 0, 1, &undirected_copy;
1971    goto &is_connected;
1972}
1973
1974*weakly_connected = \&is_weakly_connected;
1975
1976sub weakly_connected_components {
1977    &expect_directed;
1978    splice @_, 0, 1, &undirected_copy;
1979    goto &connected_components;
1980}
1981
1982sub weakly_connected_component_by_vertex {
1983    &expect_directed;
1984    splice @_, 0, 1, &undirected_copy;
1985    goto &connected_component_by_vertex;
1986}
1987
1988sub weakly_connected_component_by_index {
1989    &expect_directed;
1990    splice @_, 0, 1, &undirected_copy;
1991    goto &connected_component_by_index;
1992}
1993
1994sub same_weakly_connected_components {
1995    &expect_directed;
1996    splice @_, 0, 1, &undirected_copy;
1997    goto &same_connected_components;
1998}
1999
2000sub weakly_connected_graph {
2001    &expect_directed;
2002    splice @_, 0, 1, &undirected_copy;
2003    goto &connected_graph;
2004}
2005
2006sub _strongly_connected_components_compute {
2007    my $g = $_[0];
2008    require Graph::Traversal::DFS;
2009    require List::Util;
2010    my $t = Graph::Traversal::DFS->new($g);
2011    my @d = reverse $t->dfs;
2012    my @c;
2013    my %v2c;
2014    my $u = Graph::Traversal::DFS->new(
2015	$g->transpose_graph,
2016	next_root => sub {
2017	    my ($t, $u) = @_;
2018	    return if !defined(my $root = List::Util::first(
2019		sub { exists $u->{$_} }, @d
2020	    ));
2021	    push @c, [];
2022	    return $root;
2023	},
2024	pre => sub {
2025	    my ($v, $t) = @_;
2026	    push @{ $c[-1] }, $v;
2027	    $v2c{$v} = $#c;
2028	},
2029	next_alphabetic => 1,
2030    );
2031    $u->dfs;
2032    return [ \@c, \%v2c ];
2033}
2034
2035sub _strongly_connected_components_v2c {
2036    &_strongly_connected_components->[1];
2037}
2038
2039sub _strongly_connected_components_arrays {
2040    @{ &_strongly_connected_components->[0] };
2041}
2042
2043sub _strongly_connected_components {
2044    _check_cache($_[0], 'strong_connectivity', [],
2045			   \&_strongly_connected_components_compute);
2046}
2047
2048sub strongly_connected_components {
2049    &expect_directed;
2050    goto &_strongly_connected_components_arrays;
2051}
2052
2053sub strongly_connected_component_by_vertex {
2054    &expect_directed;
2055    &_strongly_connected_components_v2c->{$_[1]};
2056}
2057
2058sub strongly_connected_component_by_index {
2059    &expect_directed;
2060    my $i = $_[1];
2061    return if !defined(my $c = &_strongly_connected_components->[0][ $i ]);
2062    @$c;
2063}
2064
2065sub same_strongly_connected_components {
2066    &expect_directed;
2067    my ($g, @args) = @_;
2068    require Set::Object;
2069    Set::Object->new(@{ &_strongly_connected_components_v2c }{@args})->size <= 1;
2070}
2071
2072sub is_strongly_connected {
2073    &strongly_connected_components == 1;
2074}
2075
2076*strongly_connected = \&is_strongly_connected;
2077
2078sub strongly_connected_graph {
2079    &expect_directed;
2080    my ($g, %attr) = @_;
2081    my $sc_cb = \&_super_component;
2082    _opt_get(\%attr, super_component => \$sc_cb);
2083    _opt_unknown(\%attr);
2084    my ($c, $v2c) = @{ &_strongly_connected_components };
2085    my $s = Graph->new;
2086    my @s = map $sc_cb->(@$_), @$c;
2087    $s->set_vertex_attribute($s[$_], 'subvertices', $c->[$_]) for 0..$#$c;
2088    require List::Util;
2089    $s->add_edges(map [@s[ @$v2c{ @$_ } ]], grep List::Util::uniq( @$v2c{ @$_ } ) > 1, &_edges05);
2090    return $s;
2091}
2092
2093###
2094# Biconnectivity.
2095#
2096
2097sub _biconnectivity_out {
2098  my ($state, $u, $v) = @_;
2099  my @BC;
2100  while (@{$state->{stack}}) {
2101    push @BC, my $e = pop @{$state->{stack}};
2102    last if $e->[0] eq $u && $e->[1] eq $v;
2103  }
2104  push @{$state->{BC}}, \@BC if @BC;
2105}
2106
2107sub _biconnectivity_dfs {
2108  my ($g, $u, $state) = @_;
2109  $state->{low}{$u} = $state->{num}{$u} = $state->{dfs}++;
2110  for my $v ($g->successors($u)) {
2111    if (!exists $state->{num}{$v}) {
2112      push @{$state->{stack}}, [$u, $v];
2113      $state->{pred}{$v} = $u;
2114      $state->{succ}{$u}{$v}++;
2115      _biconnectivity_dfs($g, $v, $state);
2116      $state->{low}{$u} = List::Util::min(@{ $state->{low} }{$u, $v});
2117      _biconnectivity_out($state, $u, $v)
2118	if $state->{low}{$v} >= $state->{num}{$u};
2119    } elsif (defined $state->{pred}{$u} &&
2120	     $state->{pred}{$u} ne $v &&
2121	     $state->{num}{$v} < $state->{num}{$u}) {
2122      push @{$state->{stack}}, [$u, $v];
2123      $state->{low}{$u} = List::Util::min($state->{low}{$u}, $state->{num}{$v});
2124    }
2125  }
2126}
2127
2128sub _biconnectivity_compute {
2129    require List::Util;
2130    my ($g) = @_;
2131    my %state = (BC=>[], dfs=>0);
2132    my @u = $g->vertices;
2133    for my $u (@u) {
2134	next if exists $state{num}->{$u};
2135	_biconnectivity_dfs($g, $u, \%state);
2136	push @{$state{BC}}, delete $state{stack} if @{ $state{stack} || _empty_array };
2137    }
2138
2139    # Mark the components each vertex belongs to.
2140    my ($bci, %v2bc, %bc2v) = 0;
2141    for my $bc (@{$state{BC}}) {
2142      $v2bc{$_}{$bci} = undef for map @$_, @$bc;
2143      $bci++;
2144    }
2145
2146    # Any isolated vertices get each their own component.
2147    $v2bc{$_}{$bci++} = undef for grep !exists $v2bc{$_}, @u;
2148
2149    # build vector now we know how big to make it
2150    my ($Z, %v2bc_vec, @ap) = "\0" x (($bci + 7) / 8);
2151    @v2bc_vec{@u} = ($Z) x @u;
2152    for my $v (@u) {
2153	my @components = keys %{ $v2bc{$v} };
2154	vec($v2bc_vec{$v}, $_, 1) = 1 for @components;
2155	$bc2v{$_}{$v}{$_} = undef for @components;
2156	# Articulation points / cut vertices are the vertices
2157	# which belong to more than one component.
2158	push @ap, $v if @components > 1;
2159    }
2160
2161    # Bridges / cut edges are the components of two vertices.
2162    my @br = grep @$_ == 2, map [keys %$_], values %bc2v;
2163
2164    # Create the subgraph components.
2165    my @sg = map [ List::Util::uniq( map @$_, @$_ ) ], @{$state{BC}};
2166    return [ \@ap, \@sg, \@br, \%v2bc, \%v2bc_vec, $Z ];
2167}
2168
2169sub biconnectivity {
2170    &expect_undirected;
2171    @{ _check_cache($_[0], 'biconnectivity', [],
2172			   \&_biconnectivity_compute, @_[1..$#_]) || _empty_array };
2173}
2174
2175sub is_biconnected {
2176    &edges >= 2 ? @{ (&biconnectivity)[0] } == 0 : undef ;
2177}
2178
2179sub is_edge_connected {
2180    &edges >= 2 ? @{ (&biconnectivity)[2] } == 0 : undef;
2181}
2182
2183sub is_edge_separable {
2184    &edges >= 2 ? @{ (&biconnectivity)[2] } > 0 : undef;
2185}
2186
2187sub articulation_points {
2188    @{ (&biconnectivity)[0] };
2189}
2190
2191*cut_vertices = \&articulation_points;
2192
2193sub biconnected_components {
2194    @{ (&biconnectivity)[1] };
2195}
2196
2197sub biconnected_component_by_index {
2198    my ($i) = splice @_, 1, 1;
2199    (&biconnectivity)[1]->[ $i ];
2200}
2201
2202sub biconnected_component_by_vertex {
2203    my ($v) = splice @_, 1, 1;
2204    my $v2bc = (&biconnectivity)[3];
2205    splice @_, 1, 0, $v;
2206    return defined $v2bc->{ $v } ? keys %{ $v2bc->{ $v } } : ();
2207}
2208
2209sub same_biconnected_components {
2210    my ($v2bc, $Z) =  (&biconnectivity)[4,5];
2211    return 0 if grep !defined, my @vecs = @$v2bc{ @_[1..$#_] };
2212    my $accumulator = $vecs[0];
2213    $accumulator &= $_ for @vecs[1..$#vecs]; # accumulate 0s -> all in same
2214    $accumulator ne $Z;
2215}
2216
2217sub biconnected_graph {
2218    my ($g, %opt) = @_;
2219    my $bc = (&biconnectivity)[1];
2220    my $bcg = Graph->new(directed => 0);
2221    my $sc_cb = $opt{super_component} || \&_super_component;
2222    my @s = map $sc_cb->(@$_), @$bc;
2223    $bcg->set_vertex_attribute($s[$_], 'subvertices', $bc->[$_]) for 0..$#$bc;
2224    my @edges;
2225    for my $i (0..$#$bc) {
2226	my @u = @{ $bc->[ $i ] };
2227	for my $j (0..$i-1) {
2228	    my %j; @j{ @{ $bc->[ $j ] } } = ();
2229	    next if !grep exists $j{ $_ }, @u;
2230	    push @edges, [ @s[$i, $j] ];
2231	}
2232    }
2233    $bcg->add_edges(@edges);
2234    return $bcg;
2235}
2236
2237sub bridges {
2238    @{ (&biconnectivity)[2] || _empty_array };
2239}
2240
2241###
2242# SPT.
2243#
2244
2245sub _SPT_add {
2246    my ($g, $h, $HF, $r, $attr, $unseen, $etc) = @_;
2247    my $etc_r = $etc->{ $r } || 0;
2248    for my $s ( grep exists $unseen->{ $_ }, $g->successors( $r ) ) {
2249	my $t = $g->get_edge_attribute( $r, $s, $attr );
2250	$t = 1 unless defined $t;
2251	__carp_confess "Graph::SPT_Dijkstra: edge $r-$s is negative ($t)"
2252	    if $t < 0;
2253	if (!defined($etc->{ $s }) || ($etc_r + $t) < $etc->{ $s }) {
2254	    my $etc_s = $etc->{ $s } || 0;
2255	    $etc->{ $s } = $etc_r + $t;
2256	    # print "$r - $s : setting $s to $etc->{ $s } ($etc_r, $etc_s)\n";
2257	    $h->set_vertex_attributes($s, { $attr=>$etc->{ $s }, 'p', $r });
2258	    $HF->add( Graph::SPTHeapElem->new($r, $s, $etc->{ $s }) );
2259	}
2260    }
2261}
2262
2263sub _SPT_Dijkstra_compute {
2264    require Graph::SPTHeapElem;
2265    my $sptg = $_[0]->_heap_walk($_[0]->new, \&_SPT_add, {}, @_[1..$#_]);
2266    $sptg->set_graph_attribute('SPT_Dijkstra_root', $_[4]);
2267    $sptg;
2268}
2269
2270sub SPT_Dijkstra {
2271    my $g = $_[0];
2272    my @args = &_root_opt;
2273    _check_cache($g, 'SPT_Dijkstra', [$args[3]],
2274	\&_SPT_Dijkstra_compute, @args);
2275}
2276
2277*SSSP_Dijkstra = \&SPT_Dijkstra;
2278
2279*single_source_shortest_paths = \&SPT_Dijkstra;
2280
2281sub SP_Dijkstra {
2282    my ($g, $u, $v) = @_;
2283    my $sptg = $g->SPT_Dijkstra(first_root => $u);
2284    my @path = ($v);
2285    require Set::Object;
2286    my $seen = Set::Object->new;
2287    my $V = $g->vertices;
2288    my $p;
2289    while (defined($p = $sptg->get_vertex_attribute($v, 'p'))) {
2290	last if $seen->contains($p);
2291	push @path, $p;
2292	$v = $p;
2293	$seen->insert($p);
2294	last if $seen->size == $V || $u eq $v;
2295    }
2296    return if !@path or $path[-1] ne $u;
2297    return reverse @path;
2298}
2299
2300sub __SPT_Bellman_Ford {
2301    my ($g, $u, $v, $attr, $d, $p, $c0, $c1) = @_;
2302    return unless $c0->{ $u };
2303    my $w = $g->get_edge_attribute($u, $v, $attr);
2304    $w = 1 unless defined $w;
2305    if (defined $d->{ $v }) {
2306	if (defined $d->{ $u }) {
2307	    if ($d->{ $v } > $d->{ $u } + $w) {
2308		$d->{ $v } = $d->{ $u } + $w;
2309		$p->{ $v } = $u;
2310		$c1->{ $v }++;
2311	    }
2312	} # else !defined $d->{ $u } &&  defined $d->{ $v }
2313    } else {
2314	if (defined $d->{ $u }) {
2315	    #  defined $d->{ $u } && !defined $d->{ $v }
2316	    $d->{ $v } = $d->{ $u } + $w;
2317	    $p->{ $v } = $u;
2318	    $c1->{ $v }++;
2319	} # else !defined $d->{ $u } && !defined $d->{ $v }
2320    }
2321}
2322
2323sub _SPT_Bellman_Ford {
2324    my ($g, $opt, $unseenh, $unseena, $r, $next, $code, $attr) = @_;
2325    my %d;
2326    return unless defined $r;
2327    $d{ $r } = 0;
2328    my %p;
2329    my $V = $g->vertices;
2330    my %c0; # Changed during the last iteration?
2331    $c0{ $r }++;
2332    for (my $i = 0; $i < $V; $i++) {
2333	my %c1;
2334	for my $e ($g->edges) {
2335	    my ($u, $v) = @$e;
2336	    __SPT_Bellman_Ford($g, $u, $v, $attr, \%d, \%p, \%c0, \%c1);
2337	    __SPT_Bellman_Ford($g, $v, $u, $attr, \%d, \%p, \%c0, \%c1)
2338		if $g->undirected;
2339	}
2340	%c0 = %c1 unless $i == $V - 1;
2341    }
2342
2343    for my $e ($g->edges) {
2344	my ($u, $v) = @$e;
2345	if (defined $d{ $u } && defined $d{ $v }) {
2346	    my $d = $g->get_edge_attribute($u, $v, $attr);
2347	    __carp_confess "Graph::SPT_Bellman_Ford: negative cycle exists"
2348		if defined $d && $d{ $v } > $d{ $u } + $d;
2349	}
2350    }
2351
2352    return (\%p, \%d);
2353}
2354
2355sub _SPT_Bellman_Ford_compute {
2356    my ($g, @args) = @_;
2357    my ($p, $d) = $g->_SPT_Bellman_Ford(@args);
2358    my $h = $g->new;
2359    for my $v (keys %$p) {
2360	my $u = $p->{ $v };
2361	$h->set_edge_attribute( $u, $v, $args[6],
2362				$g->get_edge_attribute($u, $v, $args[6]));
2363	$h->set_vertex_attributes( $v, { $args[6], $d->{ $v }, p => $u } );
2364    }
2365    $h->set_graph_attribute('SPT_Bellman_Ford_root', $args[3]);
2366    $h;
2367}
2368
2369sub SPT_Bellman_Ford {
2370    my @args = &_root_opt;
2371    _check_cache($_[0], 'SPT_Bellman_Ford', [$args[3]],
2372	\&_SPT_Bellman_Ford_compute, @args);
2373}
2374
2375*SSSP_Bellman_Ford = \&SPT_Bellman_Ford;
2376
2377sub SP_Bellman_Ford {
2378    my ($g, $u, $v) = @_;
2379    my $sptg = $g->SPT_Bellman_Ford(first_root => $u);
2380    my @path = ($v);
2381    require Set::Object;
2382    my $seen = Set::Object->new;
2383    my $V = $g->vertices;
2384    my $p;
2385    while (defined($p = $sptg->get_vertex_attribute($v, 'p'))) {
2386	last if $seen->contains($p);
2387	push @path, $p;
2388	$v = $p;
2389	$seen->insert($p);
2390	last if $seen->size == $V;
2391    }
2392    # @path = () if @path && "$path[-1]" ne "$u";
2393    return reverse @path;
2394}
2395
2396###
2397# Transitive Closure.
2398#
2399
2400sub TransitiveClosure_Floyd_Warshall {
2401    my $self = shift;
2402    require Graph::TransitiveClosure;
2403    Graph::TransitiveClosure->new($self, @_);
2404}
2405
2406*transitive_closure = \&TransitiveClosure_Floyd_Warshall;
2407
2408sub APSP_Floyd_Warshall {
2409    my $self = shift;
2410    require Graph::TransitiveClosure;
2411    Graph::TransitiveClosure->new($self, path => 1, @_);
2412}
2413
2414*all_pairs_shortest_paths = \&APSP_Floyd_Warshall;
2415
2416sub _transitive_closure_matrix_compute {
2417    &APSP_Floyd_Warshall->transitive_closure_matrix;
2418}
2419
2420sub transitive_closure_matrix {
2421    _check_cache($_[0], 'transitive_closure_matrix', [],
2422	\&_transitive_closure_matrix_compute, @_[1..$#_]);
2423}
2424
2425sub path_length {
2426    shift->transitive_closure_matrix->path_length(@_);
2427}
2428
2429sub path_successor {
2430    shift->transitive_closure_matrix->path_successor(@_);
2431}
2432
2433sub path_vertices {
2434    shift->transitive_closure_matrix->path_vertices(@_);
2435}
2436
2437sub all_paths {
2438    shift->transitive_closure_matrix->all_paths(@_);
2439}
2440
2441sub is_reachable {
2442    shift->transitive_closure_matrix->is_reachable(@_);
2443}
2444
2445sub for_shortest_paths {
2446    my $g = shift;
2447    my $c = shift;
2448    my $t = $g->transitive_closure_matrix;
2449    my @v = $g->vertices;
2450    my $n = 0;
2451    for my $u (@v) {
2452	$c->($t, $u, $_, ++$n) for grep $t->is_reachable($u, $_), @v;
2453    }
2454    return $n;
2455}
2456
2457sub _minmax_path {
2458    my $g = shift;
2459    my $min;
2460    my $max;
2461    my $minp;
2462    my $maxp;
2463    $g->for_shortest_paths(sub {
2464			       my ($t, $u, $v, $n) = @_;
2465			       my $l = $t->path_length($u, $v);
2466			       return unless defined $l;
2467			       my $p;
2468			       if ($u ne $v && (!defined $max || $l > $max)) {
2469				   $max = $l;
2470				   $maxp = $p = [ $t->path_vertices($u, $v) ];
2471			       }
2472			       if ($u ne $v && (!defined $min || $l < $min)) {
2473				   $min = $l;
2474				   $minp = $p || [ $t->path_vertices($u, $v) ];
2475			       }
2476			   });
2477    return ($min, $max, $minp, $maxp);
2478}
2479
2480sub diameter {
2481    my $g = shift;
2482    my ($min, $max, $minp, $maxp) = $g->_minmax_path(@_);
2483    return defined $maxp ? (wantarray ? @$maxp : $max) : undef;
2484}
2485
2486*graph_diameter = \&diameter;
2487
2488sub longest_path {
2489    my ($g, $u, $v) = @_;
2490    my $t = $g->transitive_closure_matrix;
2491    if (defined $u) {
2492	return wantarray ? $t->path_vertices($u, $v) : $t->path_length($u, $v)
2493	    if defined $v;
2494	my $max;
2495	my @max;
2496	for my $v (grep $u ne $_, $g->vertices) {
2497	    my $l = $t->path_length($u, $v);
2498	    next if !(defined $l && (!defined $max || $l > $max));
2499	    $max = $l;
2500	    @max = $t->path_vertices($u, $v);
2501	}
2502	return wantarray ? @max : $max;
2503    }
2504    if (defined $v) {
2505	my $max;
2506	my @max;
2507	for my $u (grep $_ ne $v, $g->vertices) {
2508	    my $l = $t->path_length($u, $v);
2509	    next if !(defined $l && (!defined $max || $l > $max));
2510	    $max = $l;
2511	    @max = $t->path_vertices($u, $v);
2512	}
2513	return wantarray ? @max : @max - 1;
2514    }
2515    my ($min, $max, $minp, $maxp) = $g->_minmax_path(@_);
2516    return defined $maxp ? (wantarray ? @$maxp : $max) : undef;
2517}
2518
2519sub vertex_eccentricity {
2520    &expect_undirected;
2521    my ($g, $u) = @_;
2522    return Infinity() if !&is_connected;
2523    my $max;
2524    for my $v (grep $u ne $_, $g->vertices) {
2525	my $l = $g->path_length($u, $v);
2526	next if !(defined $l && (!defined $max || $l > $max));
2527	$max = $l;
2528    }
2529    return defined $max ? $max : Infinity();
2530}
2531
2532sub shortest_path {
2533    &expect_undirected;
2534    my ($g, $u, $v) = @_;
2535    my $t = $g->transitive_closure_matrix;
2536    if (defined $u) {
2537	return wantarray ? $t->path_vertices($u, $v) : $t->path_length($u, $v)
2538	    if defined $v;
2539	my $min;
2540	my @min;
2541	for my $v (grep $u ne $_, $g->vertices) {
2542	    my $l = $t->path_length($u, $v);
2543	    next if !(defined $l && (!defined $min || $l < $min));
2544	    $min = $l;
2545	    @min = $t->path_vertices($u, $v);
2546	}
2547	# print "min/1 = @min\n";
2548	return wantarray ? @min : $min;
2549    }
2550    if (defined $v) {
2551	my $min;
2552	my @min;
2553	for my $u (grep $_ ne $v, $g->vertices) {
2554	    my $l = $t->path_length($u, $v);
2555	    next if !(defined $l && (!defined $min || $l < $min));
2556	    $min = $l;
2557	    @min = $t->path_vertices($u, $v);
2558	}
2559	# print "min/2 = @min\n";
2560	return wantarray ? @min : $min;
2561    }
2562    my ($min, $max, $minp, $maxp) = $g->_minmax_path(@_);
2563    return if !defined $minp;
2564    wantarray ? @$minp : $min;
2565}
2566
2567sub radius {
2568    &expect_undirected;
2569    my $g = shift;
2570    my ($center, $radius) = (undef, Infinity());
2571    for my $v ($g->vertices) {
2572	my $x = $g->vertex_eccentricity($v);
2573	($center, $radius) = ($v, $x) if defined $x && $x < $radius;
2574    }
2575    return $radius;
2576}
2577
2578sub center_vertices {
2579    &expect_undirected;
2580    my ($g, $delta) = @_;
2581    $delta = 0 unless defined $delta;
2582    $delta = abs($delta);
2583    my @c;
2584    my $Inf = Infinity();
2585    my $r = $g->radius;
2586    if (defined $r && $r != $Inf) {
2587	for my $v ($g->vertices) {
2588	    my $e = $g->vertex_eccentricity($v);
2589	    next unless defined $e && $e != $Inf;
2590	    push @c, $v if abs($e - $r) <= $delta;
2591	}
2592    }
2593    return @c;
2594}
2595
2596*centre_vertices = \&center_vertices;
2597
2598sub average_path_length {
2599    my $g = shift;
2600    my @A = @_;
2601    my $d = 0;
2602    my $m = 0;
2603    $g->for_shortest_paths(sub {
2604        my ($t, $u, $v, $n) = @_;
2605        return unless my $l = $t->path_length($u, $v);
2606        return if defined $A[0] && $u ne $A[0];
2607        return if defined $A[1] && $v ne $A[1];
2608        $d += $l;
2609        $m++;
2610    });
2611    return $m ? $d / $m : undef;
2612}
2613
2614###
2615# Simple tests.
2616#
2617
2618sub is_multi_graph {
2619    return 0 unless &is_multiedged || &is_countedged;
2620    my $g = $_[0];
2621    my $multiedges = 0;
2622    for my $e (&_edges05) {
2623	my ($u, @v) = @$e;
2624	return 0 if grep $u eq $_, @v;
2625	$multiedges++ if $g->get_edge_count(@$e) > 1;
2626    }
2627    return $multiedges;
2628}
2629
2630sub is_simple_graph {
2631    return 1 unless &is_multiedged || &is_countedged;
2632    my $g = $_[0];
2633    return 0 if grep $g->get_edge_count(@$_) > 1, &_edges05;
2634    return 1;
2635}
2636
2637sub is_pseudo_graph {
2638    my $m = &is_countedged || &is_multiedged;
2639    my $g = $_[0];
2640    for my $e (&_edges05) {
2641	my ($u, @v) = @$e;
2642	return 1 if grep $u eq $_, @v;
2643	return 1 if $m && $g->get_edge_count($u, @v) > 1;
2644    }
2645    return 0;
2646}
2647
2648###
2649# Rough isomorphism guess.
2650#
2651
2652my %_factorial = (0 => 1, 1 => 1);
2653
2654sub __factorial {
2655    my $n = shift;
2656    for (my $i = 2; $i <= $n; $i++) {
2657	next if exists $_factorial{$i};
2658	$_factorial{$i} = $i * $_factorial{$i - 1};
2659    }
2660    $_factorial{$n};
2661}
2662
2663sub _factorial {
2664    my $n = int(shift);
2665    __carp_confess "factorial of a negative number" if $n < 0;
2666    __factorial($n) unless exists $_factorial{$n};
2667    return $_factorial{$n};
2668}
2669
2670sub could_be_isomorphic {
2671    my ($g0, $g1) = @_;
2672    return 0 unless &vertices == $g1->vertices;
2673    return 0 unless &_edges05  == $g1->_edges05;
2674    my %d0;
2675    $d0{ $g0->in_degree($_) }{ $g0->out_degree($_) }++ for &vertices;
2676    my %d1;
2677    $d1{ $g1->in_degree($_) }{ $g1->out_degree($_) }++ for $g1->vertices;
2678    return 0 unless keys %d0 == keys %d1;
2679    for my $da (keys %d0) {
2680	return 0
2681	    unless exists $d1{$da} &&
2682		   keys %{ $d0{$da} } == keys %{ $d1{$da} };
2683	return 0
2684	    if grep !(exists $d1{$da}{$_} && $d0{$da}{$_} == $d1{$da}{$_}),
2685	    keys %{ $d0{$da} };
2686    }
2687    for my $da (keys %d0) {
2688	return 0 if grep $d1{$da}{$_} != $d0{$da}{$_}, keys %{ $d0{$da} };
2689	delete $d1{$da};
2690    }
2691    return 0 unless keys %d1 == 0;
2692    my $f = 1;
2693    for my $da (keys %d0) {
2694	$f *= _factorial(abs($d0{$da}{$_})) for keys %{ $d0{$da} };
2695    }
2696    return $f;
2697}
2698
2699###
2700# Analysis functions.
2701
2702sub subgraph_by_radius {
2703    $_[0]->subgraph([ @_[1..$#_-1], &reachable_by_radius ]);
2704}
2705
2706sub clustering_coefficient {
2707    my ($g) = @_;
2708    return unless my @v = $g->vertices;
2709    require Set::Object;
2710    my %clustering;
2711
2712    my $gamma = 0;
2713
2714    for my $n (@v) {
2715	my $gamma_v = 0;
2716	my @neigh = $g->successors($n);
2717	my $c = Set::Object->new;
2718	for my $u (@neigh) {
2719	    for my $v (grep +(!$c->contains("$u-$_") && $g->has_edge($u, $_)), @neigh) {
2720		$gamma_v++;
2721		$c->insert("$u-$v");
2722		$c->insert("$v-$u");
2723	    }
2724	}
2725	if (@neigh > 1) {
2726	    $clustering{$n} = $gamma_v/(@neigh * (@neigh - 1) / 2);
2727	    $gamma += $gamma_v/(@neigh * (@neigh - 1) / 2);
2728	} else {
2729	    $clustering{$n} = 0;
2730	}
2731    }
2732
2733    $gamma /= @v;
2734
2735    return wantarray ? ($gamma, %clustering) : $gamma;
2736}
2737
2738sub betweenness {
2739    my $g = shift;
2740
2741    my @V = $g->vertices();
2742
2743    my %Cb; # C_b{w} = 0
2744
2745    @Cb{@V} = ();
2746
2747    for my $s (@V) {
2748	my @S; # stack (unshift, shift)
2749
2750	my %P; # P{w} = empty list
2751	$P{$_} = [] for @V;
2752
2753	my %sigma; # \sigma{t} = 0
2754	$sigma{$_} = 0 for @V;
2755	$sigma{$s} = 1;
2756
2757	my %d; # d{t} = -1;
2758	$d{$_} = -1 for @V;
2759	$d{$s} = 0;
2760
2761	my @Q; # queue (push, shift)
2762	push @Q, $s;
2763
2764	while (@Q) {
2765	    my $v = shift @Q;
2766	    unshift @S, $v;
2767	    for my $w ($g->successors($v)) {
2768		# w found for first time
2769		if ($d{$w} < 0) {
2770		    push @Q, $w;
2771		    $d{$w} = $d{$v} + 1;
2772		}
2773		# Shortest path to w via v
2774		if ($d{$w} == $d{$v} + 1) {
2775		    $sigma{$w} += $sigma{$v};
2776		    push @{ $P{$w} }, $v;
2777		}
2778	    }
2779	}
2780
2781	my %delta;
2782	$delta{$_} = 0 for @V;
2783
2784	while (@S) {
2785	    my $w = shift @S;
2786	    $delta{$_} += $sigma{$_}/$sigma{$w} * (1 + $delta{$w})
2787		for @{ $P{$w} };
2788	    $Cb{$w} += $delta{$w} if $w ne $s;
2789	}
2790    }
2791
2792    return %Cb;
2793}
2794
27951;
2796