#!perl

=head1 NAME

hal - finds non-isomorph quartic trees

=head1 SYNOPSIS

B<hal> I<n>

=head1 DESCRIPTION

The program finds all non-isomorph unrooted trees
in which all nodes have a degree of at most 4.
Its argument I<n> is the number of nodes.

The number of such trees is the sequence
A000602 ("http://www.research.att.com/projects/OEIS?Anum=A000602"),
which starts like this:

1, 1, 1, 1, 2, 3, 5, 9, 18, 35, 75, 159, 355, 802, 1858, 4347, 10359, 
+24894, 
60523, 148284, 366319, 910726, 2278658, 5731580, 14490245, 
36797588, 93839412, 240215803, 617105614, 1590507121, 4111846763, 
10660307791, 27711253769

=head1 OUTPUT

The trees are printed in a symbolic form.  
For some reason, all nodes are represented by a letter C 
(node = csúcs in Hungarian).  Edges are represented by a minus sign.
If a node has multiple children, then all but the last are parenthisiz
+ed.

For example, this formula would correspond to this tree:

 A-B-C(-D)(-E)-F-G(-H(-I)-J)-K-L-M

     D
     |
 A-B-C-F-G-K-L-M
     |   |
     E   H-I
         |
         J

=head1 BUGS

The implementation is somewhat slow, so it works for
small integer I<n> only.

=cut

# head --------------------------------------------------------

use warnings; use strict;

use List::Util qw(min max first);
use Carp qw(confess cluck);
use Data::Dumper;


# helpers -------------------------------------------------

sub INF() { exp(exp(42)); }

sub insert1 {
# sorts a list if all but the first one is already sorted in _decreasi
+ng_ order
# returns the insert position and the new list
 my($e, @l) = @_;
 my @r = sort({ $b <=> $a } $e, @l);
 my $p = first(sub { $r[$_] == $e }, 0 .. @r - 1);
 $p, \@r;
}


# data ---------------------------------------------------

# format of an entry of @hal:
# { c => [ indices of children in @hal ], 
#   w => (number of carbon atoms, excluding the root one),
#   n => (normalized form) }

my @hal = {"c" => [], "w" => 0};

my @w_inv = 0;

my %hal_inv = (); # 0 nem jöhet ki!


# search --------------------------------------------------------

sub hal_search {
 my($n) = @_;
 for my $w (1..$n) { 
  $w_inv[$w] = @hal;
  hal_search1($w); 
  $w_inv[$w + 1] = @hal;
 }
}

sub hal_search1 {
 my($w) = @_;
 for my $w1 (0 .. $w - 1) {
  for my $w2 (0 .. $w - 1) {
   my $w3 = $w - 1 - $w1 - $w2;
   if ($w - 1 >= $w1 && $w1 >= $w2 && $w2 >= $w3 && $w3 >= 0)
    { hal_search2($w, $w1, $w2, $w3); }
  }
 }
}

sub hal_search2 {
 my($w, $w1, $w2, $w3) = @_;
 hal_loopw($w1, -1, sub {
  my($c1) = @_;
  hal_loopw($w2, $c1, sub {
   my($c2) = @_;
   hal_loopw($w3, $c2, sub {
    my($c3) = @_;
    hal_find($w, $c1, $c2, $c3);
   });
  });
 });
}

sub hal_loopw {
 my($w, $m, $f) = @_;
 my($a, $b);
 $a = $w_inv[$w];
 $b = $w_inv[$w+1] - 1;
 0 <= $m and $m < $b and $b = $m;
 for my $k ($a .. $b) { 
  &$f($k); 
 }
}

sub hal_find {
 my($w, $c1, $c2, $c3) = @_;
 $w == 1 + ${$hal[$c1]}{"w"} + ${$hal[$c2]}{"w"} + ${$hal[$c3]}{"w"} o
+r
  confess "internal error";
 push @hal, {"c" => [grep({ 0 != $_ } $c1, $c2, $c3)], "w" => $w};
 $hal_inv{join(",", grep({$_ != 0} $c1, $c2, $c3))} = @hal - 1;
 return;
}


# norm -----------------------------------------------------

sub normalize {
 my($n) = @_;
 hal_loopw($n, -1, sub {
  norm_try ($_[0]); 
 });
}

sub norm_try {
 my($n) = @_;
 my $h = $hal[$n];
 exists $$h{"n"} and return;
 norm_mark($n, -1, $n);
}

sub norm_mark {
 my($n, $xl, @c) = @_;
 @c == 1 and do {
  ${$hal[$c[0]]}{"n"} = $n;
 };
 for my $k (0 .. @c - 1) {
  $k == $xl and next;
  my $j = join(",", @c[0 .. $k - 1, $k + 1 .. @c - 1]);
  exists($hal_inv{$j}) or
   confess "internal error: can not rotate tree ($j)";
  my $a = $hal_inv{$j};
  my @n = @{${$hal[$c[$k]]}{"c"}};
  my ($p, $l) = insert1($a, @n);
  norm_mark($n, $p, @$l);
 }
}


# format --------------------------------------------------------

sub hal_strmol {
 "C" . hal_str($_[0]);
}

sub hal_str {
 my($n) = @_;
 $n == 0 and confess "internal error";
 my($c1, @c2) = @{${$hal[$n]}{"c"}};
 my $s = "-C";
 for my $k (reverse @c2) {
  $s .= hal_strp($k);
 }
 defined($c1) && $c1 != 0 and
  $s .= hal_str($c1);
 $s;
}

sub hal_strp {
 0 == $_[0] ? "" : "(" . hal_str($_[0]) . ")";
 # "(" . hal_str($_[0]) . ")";
}

sub printresults {
 my($n) = @_;
 my %b;
 hal_loopw($n, -1, sub {
  my($k) = @_;
  my $l = ${$hal[$k]}{"n"};
  $b{$l} = $k;
 });
 for my $l (sort(keys(%b))) {
  my $k = $b{$l};
  print hal_strmol($k), "\n"; 
 }
 print 0+keys(%b), "\n";
}


# main --------------------------------------------------

my $MAXN = 5 - 1;
@ARGV > 0 and $ARGV[0]=~/^(\d+)$/ and $1 >= 2 and 
 $MAXN = $1 - 1;

warn "Doing search, " . ($MAXN+1) . "...\n";
hal_search($MAXN);

warn "Finding isomorphic trees...";
normalize($MAXN);

warn "Printing results...";
printresults($MAXN);




__END__
[download]