#!/usr/bin/env perl

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

use List::Util qw(product minstr uniqstr);
use Array::Set qw(set_diff set_union set_intersect set_symdiff);
use Algorithm::Combinatorics qw(subsets permutations);

use Getopt::Long;

GetOptions(
    '--count-only'  => \my $count_only,
    '--no-symmetry' => \my $no_symmetry,
    '--positive-singletons' => \my $positive_singletons,
) or die 'failed parsing options';

my $cube = Cube(shift // die 'need dimension');
my $N = $cube->set;

my @group;
for my $p (permutations $N) {
    for my $Z (subsets $N) {
        push @group, 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 to_str {
    join '', map { $_ > 0 ? '+' : '-' } shift->@*
}

my $sat = CInet::ManySAT->new;
# The empty set is always positive.
$sat->add([ $cube->pack([ [], [  ] ]) ]);
# Diamond property:
for my $A (subsets($N)) {
    next if @$A > @$N - 2;
    for my $ij (subsets(set_diff($N, $A), 2)) {
        # CNF of (A = ijA) => (iA = jA) for boolean variables ijA, iA, jA, A.
        my $a   = $cube->pack([ [], $A ]);
        my $ia  = $cube->pack([ [], set_union($A, [$ij->[0]]) ]);
        my $ja  = $cube->pack([ [], set_union($A, [$ij->[1]]) ]);
        my $ija = $cube->pack([ [], set_union($A, [$ij->@* ]) ]);
        $sat->add([  $a,  $ija,  $ia, -$ja ]);
        $sat->add([  $a,  $ija, -$ia,  $ja ]);
        $sat->add([ -$a, -$ija,  $ia, -$ja ]);
        $sat->add([ -$a, -$ija, -$ia,  $ja ]);
    }
}

# Up to the hyperoctahedral symmetry, we can always assume that the
# singletons are positive.
if ($positive_singletons) {
    $sat->add([ $cube->pack([ [], [ $_ ] ]) ]) for @$N;
}

say STDERR $sat->count(risk => 0.01), ' compatible sign vectors';
exit if $count_only;

my $all = $sat->all;
if ($no_symmetry) {
    while (defined(my $v = $all->next)) {
        say "$v";
    }
    exit;
}

print STDERR 'Performing symmetry reduction ';
my (%orbits, %reps);
while (defined(my $v = $all->next)) {
    my $V = [ map { $_ > 0 ? 1 : -1 } $v->@* ];
    next if defined $reps{to_str $V};
    my @orbit = uniqstr sort map { to_str act($V, $_) } @group;
    my $rep = minstr @orbit;
    $reps{$_} = $rep for @orbit;
    push $orbits{$rep}->@*, $_ for @orbit;
    print STDERR '.';
}
print STDERR "\n";

say STDERR scalar keys %orbits, ' representatives';
say "$_: @{[ 0+ $orbits{$_}->@* ]}" for sort keys %orbits;
