#!/usr/bin/env perl

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

use Algorithm::Combinatorics qw(subsets);
use Math::Polynomial::Multivariate;
use Matrix::Simple;

GetOptions; # process -- for argument separation

my $cube = Cube([0..4]);
my $pattern = shift // die 'need pattern';

sub c { Math::Polynomial::Multivariate->const(shift) }
sub v { Math::Polynomial::Multivariate->var(shift)   }

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

sub pr {
    my ($m, $K) = @_;
    det submatrix($m, $K, $K)
}

sub fmt_poly {
    my $p = shift;
    my $s = '';
    for ($p->flatten->as_monomials) {
        my ($c, $v) = @$_;
        $s .= sprintf("%+d ", $c);
        $s .= $v->{$_} > 1 ? "$_^$v->{$_} " : "$_ " for sort keys %$v;
    }
    $s =~ s/\s+$//r;
}

my $m = [
    [ c(1) , v('a'), v('b'), v('c'), v('d')],
    [v('a'),  c(1) , v('e'), v('f'), v('g')],
    [v('b'), v('e'),  c(1) , v('h'), v('i')],
    [v('c'), v('f'), v('h'),  c(1) , v('j')],
    [v('d'), v('g'), v('i'), v('j'),  c(1) ]
];

for my $i ($cube->set->@*) {
    my $j = $cube->pack([ [], [$i] ]) - 1;
    $m->[$i][$i] = c(substr($pattern,$j,1) . '1');
}

say 'Maximize';
say 'obj: +1';
say 'Subject to';
for my $l (2 .. $cube->dim) {
    for my $K (subsets($cube->set, $l)) {
        my $k = $cube->pack([ [], $K ]) - 1;
        my $p = pr($m, $K, $K);
        my $s = substr($pattern,$k,1);
        next if $s eq '*';
        $p = -$p if $s eq '-';
        say sprintf('  p%s: %s > 0.01', join('', map { $_+1 } @$K), fmt_poly $p);
    }
}
say 'Bounds';
for my $v (pr($m, $cube->set)->variables) {
    say sprintf('-200 <= %s <= 200', $v);
}
say 'End';
