#!/usr/bin/env perl

# This program takes as input the output of sample.sage.
# It accumulates all records into hyperoctahedral orbits and
# prints one record for each canonical representative of the
# orbit.

use Modern::Perl 2018;
use CInet::Base;

use Math::BigRat lib => 'GMP';
use List::Util qw(product minstr uniqstr);
use List::MoreUtils qw(firstidx);
use Array::Set qw(set_diff set_union set_intersect set_symdiff);
use Algorithm::Combinatorics qw(subsets permutations);

my $n = shift // die 'need dimension';
my $cube = Cube([ 0 .. $n-1 ]);
my $N = $cube->set;

my (@Bn, @Sn);
for my $p (permutations $N) {
    for my $Z (subsets $N) {
        push @Sn, [$p, 0+ @Bn] if @$Z == 0;
        push @Bn, lift($cube, $p, $Z);
    }
}

# Preprocess an element of the hyperoctahedral group to speed up
# the action on sign patterns.
sub lift {
    my ($cube, $p, $Z) = @_;
    my @Pn = subsets($cube->set);
    my $lifted = [ map 0, @Pn ];

    my $z = $cube->pack([ [], $Z ]) - 1;
    for my $K (@Pn) {
        my $k  = $cube->pack( $cube->permute($p => [ [], $K ]) ) - 1;
        my $zk = $cube->pack([ [], set_symdiff($Z, $K) ]) - 1;
        my $c = (set_intersect($Z, $K)->@* % 2) == 0 ? +1 : -1;
        # The new value at $k is the product of the old values at $z and $zk
        # and the sign $c.
        $lifted->[$k] = [ $z, $zk, $c ];
    }
    $lifted
}

sub act {
    my ($V, $lifted) = @_;
    [ map { $_->[2] * product($V->@[$_->@[0, 1]]) } @$lifted ]
}

sub matrix {
    my $n = shift;
    my $A = [ map [ map 0, 1 .. $n ], 1 .. $n ];
    for (my $i = 0; $i < $n; $i++) {
        for (my $j = $i; $j < $n; $j++) {
            $A->[$i][$j] = $A->[$j][$i] = Math::BigRat->new(shift);
        }
    }
    $A
}

sub submatrix {
    my ($A, $I, $J) = @_;
    [ map { [ $_->@[@$J] ] } $A->@[@$I] ]
}

sub entries {
    my $A = shift;
    my $n = @$A;
    map { my $i = $_; map $A->[$i][$_], $i .. $n-1 } 0 .. $n-1;
}

sub to_str {
    join '', map { $_ > 0 ? '+' : '-' } shift->@*
}

my %found;
my %reps;
while (<<>>) {
    chomp;
    die "can't parse record '$_'" unless /(.+): (\d+) \[(.+)\]/;
    my ($s, $count, $nums) = ($1, $2, $3);
    my $V = [map $_ . '1', split //, $s];
    my $rep = do {
        if (not exists $reps{$s}) {
            my @orbit = uniqstr sort map { to_str act($V, $_) } @Bn;
            my $rep = minstr @orbit;
            $reps{$_} = $rep for @orbit;
        }
        $reps{$s}
    };

    if (not exists $found{$rep}) {
        $found{$rep} = [ $count, undef ];
    }
    else {
        $found{$rep}->[0] += $count;
    }

    # Try to find a representation of $rep from the matrices we have
    # by a permutation. This avoids implementing the whole hyperoctahedral
    # action on matrices and is good enough with enough samples.
    if (not defined $found{$rep}->[1]) {
        for (@Sn) {
            my ($p, $idx) = @$_;
            my $lifted = $Bn[$idx];
            if ($rep eq to_str act($V, $lifted)) {
                my $A = matrix($n, split /, /, $nums);
                my @pinv = map { my $i = $_; firstidx { $i == $_ } @$p } $cube->set->@*;
                my $Ap = submatrix($A, \@pinv, \@pinv);
                $found{$rep}->[1] = join ', ', entries($Ap);
                last;
            }
        }
    }
}

for my $k (sort keys %found) {
    if (not defined $found{$k}->[1]) {
        say STDERR "No representation available for $k. Omitting...";
        next;
    }
    say "$k: @{[ $found{$k}->[0] ]} [@{[ $found{$k}->[1] ]}]";
}