Graph::ModularDecomposition - Modular decomposition of directed graphs

    use Graph::ModularDecomposition qw(pairstring_to_graph tree_to_string);
    my $g = new Graph::ModularDecomposition;

    my $h = $g->pairstring_to_graph( 'ab,ac,bc' );
    print "yes\n" if check_transitive( $h );
    print "yes\n" if $h->check_transitive; # same thing
    my $m = $h->modular_decomposition_EGMS;
    print tree_to_string( $m );

package Graph::ModularDecomposition;

use 5.006;
use strict;
use warnings;

require Exporter;
our $VERSION = '0.13';

use Graph 0.20105;
require Graph::Directed;

# NB! Exporter must come before Graph::Directed in @ISA
our @ISA = qw(Exporter Graph::Directed);

# This allows declaration	use Graph::ModularDecomposition ':all';
# may want tree_to_string, should move into own Tree::... module some day
# other exports are most likely for internal use only
# all other functions should be accessed as methods
our %EXPORT_TAGS = ( 'all' => [ qw(
	setminus
	setunion
	pairstring_to_graph
	partition_to_string
	tree_to_string
) ] );

our @EXPORT_OK = ( @{ $EXPORT_TAGS{'all'} } );

our @EXPORT = qw(
);

use Carp;

my $VSEP = '|'; # string used to separate vertices
my $WSEP = '\|'; # regexp used to separate vertices
my $PSEP = '\+'; # regexp used to separate elements of partition
my $QSEP = '+'; # string used to separate elements of partition

my $Debug = 0;

sub debug {
    my $class = shift;
    if ( ref($class) ) { $class = ref($class) }
    $Debug = shift;
    carp 'Turning ', ($Debug ? 'on' : 'off'), ' ',
	$class, ' debugging', ($Debug ? ", level $Debug" : '');
}

my $Canonical_form = 0;

sub canonical_form {
    my $class = shift;
    if ( ref($class) ) { $class = ref($class) }
    my $cf = shift;
    return $Canonical_form unless defined $cf;
    $Canonical_form = $cf;
}

sub new {
    my $self = shift;
    my $class = ref($self) ? ref($self) : $self;
    return bless $class->SUPER::new(@_), $class;
}

sub pairstring_to_graph {
    my $class = shift;
    if ( ref($class) ) { $class = ref($class) }
    my $pairs = shift;
    my $g = new $class;
    my ($p, $q);
    my $s = ( ( index( $pairs, '-' ) >= 0 ) ? '\-' : '' );
    foreach my $r ( split /,\s*/, $pairs ) {
	( $p, $q ) = split $s, $r;
	print "p=$p, q=$q\n" if $Debug > 2;
	if ( $q ) {
	    $g = $g->add_edge( $p, $q ) unless $g->has_edge( $p, $q );
	} else {
	    $g = $g->add_vertex( $p ) unless $g->has_vertex( $p );
	}
    }
    return bless $g, $class;
}

sub check_transitive {
    my $g = shift;
    my $g2 = $g->copy;
    my $h = $g->TransitiveClosure_Floyd_Warshall;
    # get rid of loops
    foreach ( $h->vertices ) { $h->delete_edge( $_, $_ ) }
    foreach ( $g2->vertices ) { $g2->delete_edge( $_, $_ ) }
    print STDERR "gdct: ", $g, ' vs. ', $h, "\n" if $Debug;
    return $h eq $g2;
}

sub setminus {
    my $X = shift;
    my $Y = shift;
    my @X = @{$X};
    print STDERR 'setminus# ', @X, ' - ', @{$Y}, ' = ' if $Debug > 1;
    foreach my $x ( @{$Y} ) {
	@X = grep $x ne $_, @X;
    }
    print STDERR @X, "\n" if $Debug > 1;
    return @X;
}

sub setunion {
    my $X = shift;
    my $Y = shift;
    my @X = @{$X};
    print STDERR 'setunion# ', @X, ' U ', @{$Y}, ' = ' if $Debug > 1;
    foreach my $x ( @{$Y} ) {
	push @X, $x unless grep $x eq $_, @X;
    }
    print STDERR @X, "\n" if $Debug > 1;
    return sort @X;
}

sub restriction {
    my $G = shift;
    if ( $Debug > 2 ) { print STDERR 'restriction(', ref($G), ")\n" }
    my $h = ($G->copy)->delete_vertices( setminus( [$G->vertices], [@_] ) );
    if ( $Debug > 1 ) {
	print STDERR 'restriction(', $G, '|', join($QSEP, @_), ') = ', $h, "\n"
    }
    return $h;
}

sub factor {
    my $G = shift;
    my $P = shift;
    my $GP = $G->copy;
    my $p;
    foreach my $X ( @{$P} ) {
	print STDERR "factor# X = $X\n" if $Debug > 1;
	print STDERR "factor# \@X = @$X\n" if $Debug > 1;
	my $newnode = join $VSEP, @{$X}; # turn nodes a, b, c into new node abc
	print STDERR "factor# newnode = $newnode\n" if $Debug > 1;
	my $a = ${$X}[0];
	print STDERR "factor# representative node $a\n" if $Debug > 1;
	if ( $newnode ne $a ) { # do nothing if singleton
	    $GP->add_vertex( $newnode );
	    foreach $p ( $GP->predecessors( $a ) ) {
		print STDERR "factor# predecessor $p\n" if $Debug > 2;
		$GP = $GP->add_edge( $p, $newnode )
		    unless $GP->has_edge( $p, $newnode );
	    }
	    foreach $p ( $GP->successors( $a ) ) {
		print STDERR "factor# successor $p\n" if $Debug > 2;
		$GP = $GP->add_edge( $newnode, $p )
		    unless $GP->has_edge( $newnode, $p );
	    }
	    $GP = $GP->delete_vertices( @{$X} );
	}
    }
    return $GP;
}

sub partition_subsets {
    my $G = shift;
    my $S = shift;
    my $w = shift;

    print STDERR 'p..n_subsets# @S = ', @{$S}, ", w = $w \n" if $Debug > 1;
    my (@A, @B, @C, @D);
    foreach my $x ( @{$S} ) {
	print STDERR 'p..n_subsets# xw = ', $x, $w if $Debug > 2;
	if ( $G->has_edge( $w, $x ) ) {
	    if ( $G->has_edge( $x, $w ) ) { # xw wx (not poset)
		push @A, $x;
		print STDERR ' A = ', @A, "\n" if $Debug > 2;
	    } else { # ~xw wx
		push @B, $x;
		print STDERR ' B = ', @B, "\n" if $Debug > 2;
	    }
	} else {
	    if ( $G->has_edge( $x, $w ) ) { # xw ~wx
		push @C, $x;
		print STDERR ' C = ', @C, "\n" if $Debug > 2;
	    } else { # ~xw ~wx
		push @D, $x;
		print STDERR ' D = ', @D, "\n" if $Debug > 2;
	    }
	}
    }
    return grep @{$_}, (\@A, \@B, \@C, \@D);
}

sub partition {
    my $G = shift;
    my $v = shift;

    print STDERR 'partition# G = ', $G, ", v = $v\n" if $Debug > 1;
    my (%L, @done, $tempset, $S, @ZS, $w);
    $S = [ setminus( [ $G->vertices ], [ $v ] ) ];
    print STDERR 'partition# @S = ', @{$S}, "\n" if $Debug > 1;
    $L{$S} = [ $v ];
    my @todo = ( $S );
    print STDERR 'partition# L{S}[0] = ', $L{$S}[0], "\n" if $Debug > 1;
    while ( @todo ) {
	$S = shift @todo;
	@ZS = @{$L{$S}};
	$w = $ZS[0];
	print STDERR 'partition# ZS = ', @ZS, "\n" if $Debug > 1;
	delete $L{$S};
	foreach my $W ( $G->partition_subsets( $S, $w ) ) {
	    print STDERR 'partition# W = ', @{$W}, "\n" if $Debug > 1;
	    $tempset = [ setunion( [ setminus( $S, $W ) ],
				[ setminus( \@ZS, [ $w ] ) ] ) ];
	    if ( @{$tempset} ) {
		print STDERR 'partition# tempset = ', @{$tempset}, "\n"
		    if $Debug > 1;
		$L{$W} = $tempset;
		push @todo, $W;
	    } else {
		push @done, $W;
	    }
	}
    }
    return \@done;
}

sub distinguishes {
    my ($g,$x,$y,$z) = @_;
    print STDERR " $x$y?", $g->has_edge($x,$y) if $Debug > 1;
    print STDERR " $x$z?", $g->has_edge($x,$z) if $Debug > 1;
    print STDERR " $y$x?", $g->has_edge($y,$x) if $Debug > 1;
    print STDERR " $z$x?", $g->has_edge($z,$x) if $Debug > 1;
    my $ret =  ( $g->has_edge($x,$y) != $g->has_edge($x,$z) )
	    || ( $g->has_edge($y,$x) != $g->has_edge($z,$x) );
    print STDERR "=$ret\n" if $Debug > 1;
    return $ret;
}

sub G {
    my $g = shift;
    my $v = shift;
    my $G = new ref($g);
    print STDERR 'G([', $g, "], $v) =...\n" if $Debug;
X:  foreach my $x ( $g->vertices ) {
	next X if ( $v eq $x );
	print STDERR 'X=', $x, "\n" if $Debug > 1;
	$G = $G->add_vertex( $x );
Y:	foreach my $y ( $g->vertices ) {
	    next Y if ( $v eq $y or $x eq $y );
	    print STDERR 'Y=', $y, "\n" if $Debug > 1;
	    if ( $g->distinguishes( $x, $y, $v ) ) {
		$G = $G->add_edge( $x, $y ) unless $G->has_edge( $x, $y );
	    }
	}
    }
    print STDERR '...G()=', $G, "\n" if $Debug;
    return $G;
}

sub tree_to_string {
    my $t = shift;
    my $s = '';
    return $s unless defined $t->{type};
    $s .= $t->{type} if $t->{type} ne 'leaf';
    $s .= '_' . $t->{col} if ( $t->{type} eq 'complete' );
    $s .= '[' . $t->{value} . ']';
    if ( $t->{type} ne 'leaf' ) {
	my $sep = '';
	$s .= '(';
	foreach ( @{$t->{children}} ) {
	    $s .= $sep . tree_to_string( $_ );
	    $sep = ';';
	}
	$s .= ')';
    }
    return $s;
}

sub partition_to_string {
    return join ',', sort (map { join $QSEP, sort @{$_} } @{+shift});
}

sub modular_decomposition_EGMS {
    my $g = shift;
    my $md = 0;
    $md ++;
    my $B = ' 'x$md;
    print STDERR $B, 'MD(', $g, ")=...\n" if $Debug;
    my $v = ($g->vertices)[0];
    print STDERR $B, 'v=', (defined($v) ? $v : 'undef'), "\n" if $Debug;

    my $t = {};
    unless ( $v ) {
	print STDERR $B, '...MD=', tree_to_string( $t ), "\n" if $Debug;
	$md --;
	return $t;
    }
    $t->{type} = 'leaf';
    $t->{children} = [];
    if ($g->canonical_form()) {
	$t->{value} = join($VSEP, sort($g->vertices));
    } else {
	$t->{value} = join($VSEP, $g->vertices);
    }
    $t->{col} = '0';

    if ( scalar $g->vertices == 1 ) {
	print STDERR $B, '...MD=', tree_to_string( $t ), "\n" if $Debug;
	$md --;
	return $t;
    }

    my $p = partition( $g, $v );
    push @{$p}, [ $v ];
    my $gd = $g->factor( $p );
    print STDERR $B, 'gd = ', $gd, "\n" if $Debug;
    my $Gdd = $gd->G($v)->strongly_connected_graph;
    print STDERR $B, 'Gdd = [', $Gdd, '], ', scalar $Gdd->vertices, "\n" if $Debug;

    my $u = $t;
    my @f;
    while ( @f = grep( $Gdd->out_degree($_) == 0 , $Gdd->vertices ) ) {
	print STDERR $B, "\@f=[@f]\n" if $Debug;
	my @s;
	foreach my $s ( $Gdd->vertices ) {
	    push @s, split(/$PSEP/, $s);
	}
	if ($g->canonical_form()) {
	    $u->{value} = join('', sort($v, @s));
	} else {
	    $u->{value} = join('', ($v, @s));
	}
	my $w = {};
	$w->{type} = 'leaf';
	$w->{children} = [];
	$w->{value} = $v;
	$w->{col} = '0';
	push @{$u->{children}}, $w;

	$Gdd->delete_vertices( @f );
	my @F;
	foreach my $f ( @f ) {
	    foreach my $F ( split /$PSEP/, $f ) {
		push @F, $F unless grep $F eq $_, @F;
	    }
	}
	print STDERR $B, "\@F=@F\n" if $Debug;
	if ( @f == 1 and @F > 1 ) {
	    $u->{type} = 'primitive';
	    $u->{col} = '0';
	} else {
	    my $x = substr $F[0], 0, 1; # single-char vertex names!
	    if ( $g->has_edge($v, $x) == $g->has_edge($x, $v) ) {
		$u->{type} = 'complete'; # 0 parallel, 1 series
		$u->{col} = $g->has_edge($v, $x) ? '1' : '0';
	    } else {
		$u->{type} = 'linear';
		$u->{col} = '01';
	    }
	}
	print STDERR $B, 'u = ', tree_to_string( $u ), "\n" if $Debug;
	foreach my $X ( @F ) {
	    my $m = $g->restriction( split /$WSEP/, $X )
		    ->modular_decomposition_EGMS;
	    if ( defined $m->{col}
		and ( $u->{col} eq $m->{col} )
		and (
		    ( $u->{type} eq 'complete' and $m->{type} eq 'complete' )
		   or ( $u->{type} eq 'linear' and $m->{type} eq 'linear' )
		)
	    ) {
		if ( $Debug ) {
		    print STDERR $B, "u->children= @{$u->{children}}\n";
		    print STDERR $B, 'm->children= ';
		    my $sep = '';
		    foreach ( @{$m->{children}} ) {
			print STDERR $sep, '[', tree_to_string( $_ ), ']';
			$sep = ', ';
		    }
		    print STDERR "\n";
		}
		push @{$u->{children}}, @{$m->{children}};
	    } else {
		push @{$u->{children}}, $m;
	    }
	}
	$u = $w;
    }
    print STDERR $B, '...MD=', tree_to_string( $t ), "\n" if $Debug;
    $md --;
    return $t;
}

sub classify {
    my $g = shift;
    return 'n' unless $g->check_transitive;
    my $s = tree_to_string( $g->modular_decomposition_EGMS );
    return 'i' if $s =~ m/^primitive\[[^\]]+\]\([^\(]*$/;
    return 'd' if $s =~ m/primitive/ and $s =~ m/complete_|linear/;
    return 's' if $s !~ m/primitive|complete_1/; # matches empty string
    return 'p' if $s =~ m/primitive|complete_1/;
    return 'u';
}

eval {require Graph::Bitvector2; 1} and # alas, circular dependency here
eval q{
        sub to_bitvector2 {
		my $g = shift;
		my @v = sort $g->vertices;
		my @bits;
		while ( @v ) {
	    	    my $x = shift @v;
	    	    foreach my $y ( @v ) {
				push @bits, (
		    		    $g->has_edge( $x, $y )
		    		    ? 1
		    		    : ( $g->has_edge( $y, $x ) ?  2 : 0 )
				);
	    	    }
		}
		return new Graph::Bitvector2 (join '', @bits);
        }
};

1;
__END__