#!/usr/bin/env perl

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

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

my $cube = Cube(4);
my $N = $cube->set;

my @group;
for my $Z (subsets($N)) {
    for my $p (permutations $N) {
        push @group, [$p, $Z];
    }
}

sub act {
    my ($V, $p, $Z) = @_;

    my $z = $cube->pack([ [], $Z ]) - 1;
    my $W = [ @$V ];
    for my $K (subsets($N)) {
        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;
        $W->[$k] = $V->[$z] * $V->[$zk] * $c;
    }
    $W
}

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 ]);
    }
}

# Compute all hyperoctahedral orbits.
my $all = $sat->all;
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;
}

# We know that all sign patterns on n=4 are representable, so we may use
# the theorem bounding the numbers of connected components recursively
# by their minors. We look at each of its eight 3-minors and take the
# maximum of their numbers of connected components. The below table is
# computed by `dimh0-pr3.pl`.

my %dimh0pr3 = (
    '+++++++-' => 4,
    '++++---+' => 4,
    '+++-----' => 4,
    '++-+----' => 4,
    '+-++----' => 4,
    '+---+--+' => 4,
    '+----+-+' => 4,
    '+-----++' => 4,
    '+++++---' => 4,
    '++++-+--' => 4,
    '++++--+-' => 4,
    '++-----+' => 4,
    '+-+----+' => 4,
    '+--+---+' => 4,
    '+---++++' => 4,
    '+-------' => 4,
    '++++++++' => 1,
    '+++-+---' => 1,
    '++-+-+--' => 1,
    '++----++' => 1,
    '+-++--+-' => 1,
    '+-+--+-+' => 1,
    '+--++--+' => 1,
    '+---+++-' => 1,
    '++++++--' => 2,
    '+++++-+-' => 2,
    '++++-++-' => 2,
    '+++----+' => 2,
    '++-+---+' => 2,
    '++------' => 2,
    '+-++---+' => 2,
    '+-+-----' => 2,
    '+--+----' => 2,
    '+---++-+' => 2,
    '+---+-++' => 2,
    '+----+++' => 2,
    '++++----' => 16,
    '+------+' => 16,
);

sub dual {
    my ($s, $cube) = @_;
    my @t;
    my $v = substr($s, $cube->pack([ [], $cube->set ])-1, 1);
    for ($cube->vertices) {
        my $c = substr($s, $cube->pack([ [], set_diff($cube->set, $_->[1]) ])-1, 1);
        my $w = "${v}1" * "${c}1";
        push @t, ($w > 0 ? '+' : '-');
    }
    join('', @t)
}

sub deletion {
    my ($s, $I) = @_;
    my @t;
    for (Cube(set_diff($N, $I))->vertices) {
        push @t, substr($s, $cube->pack($_)-1, 1);
    }
    join('', @t)
}

sub contraction {
    my ($s, $I) = @_;
    dual(deletion(dual($s, $cube), $I), Cube(set_diff($N, $I)))
}

sub ncomps {
    my $s = shift;
    return $dimh0pr3{$s} // die "sign pattern '$s' not found"
}

my $total = 0;
for my $rep (keys %orbits) {
    my $n1 = max map { ncomps deletion($rep, [$_]) } @$N;
    my $n2 = max map { ncomps contraction($rep, [$_]) } @$N;
    my $bound = max($n1, $n2);
    my $size = 0+ $orbits{$rep}->@*;
    say "$rep: at least $bound x $size = @{[ $bound * $size ]}";
    $total += $bound * $size;
}
say "Total = $total";
