Here's my idea how to do that:
#!/usr/bin/perl
use warnings;
use strict;
use Math::Polynomial::Solve qw(quadratic_roots);
use Math::Numbers;
use Data::Dumper;
my $level = 64;
my @powers_of_two = (1);
push (@powers_of_two, $powers_of_two[-1]*2) while $powers_of_two[-1] <
+ $level;
my @numbers = @ARGV;
@numbers = (1..$level-1) unless @numbers;
for my $number (@numbers) {
print "E($number/$level)";
for (1) {
my @divisors = Math::Numbers->new($number)->get_diviso
+rs();
if (grep({$_ == 2} @divisors)) {
# trivial
print " = E(" . $number/2 . "/" . $level/2 . "
+)";
}
if ($number < $level/2) {
# trivial
print " = E($number/" . $level/2 . ") & R";
}
if (@divisors > 2) {
# not a prime, decompose
my ($numerator1) = grep({$_ != 1 && $_ != $num
+ber} @divisors);
my $numerator2 = $number/$numerator1;
my $denominator1 = 2;
$denominator1 *= 2 while ($denominator1 < $num
+erator1);
my $denominator2 = $level/$denominator1;
if ($denominator2 > $numerator2) {
print " = E($numerator1/$denominator1)
+ & E($numerator2/$denominator2)";
}
}
# prime or failed to decompose nicely, try possible ov
+erlaps
for (my $overlap = 1; $overlap + $number < 2 * $level;
+ $overlap += 2) {
my @o_divisors = sort {$a <=> $b} Math::Number
+s->new($overlap)->get_divisors();
for my $numerator1 (@o_divisors) {
last if $numerator1 > $overlap/2;
my $numerator2 = $overlap/$numerator1;
# try to find integer d1, d2 such that
+ n1*n2 = level and n1*d2+n2*d1 = number + overlap
my $quad_a = $numerator1;
my $quad_b = - ($number + $overlap);
my $quad_c = $level * $numerator2;
my @roots = quadratic_roots($quad_a, $
+quad_b, $quad_c);
#print "overlap $overlap, n1 $numerato
+r1, roots " . join(',', @roots) . "\n";
for my $root (@roots) {
next if ref $root; # complex r
+oot
$root += 0;
next unless grep({$root == $_}
+ @powers_of_two);
my ($denominator1, $denominato
+r2) = sort {$a <=> $b} ($root, $level/$root);
next if $denominator1 < $numer
+ator1;
next if $denominator2 < $numer
+ator2;
print " = E($numerator1/$denom
+inator1) | E($numerator2/$denominator2)";
}
}
}
}
print "\n";
}
For a given power of two ($level = 2^n), it iterates through all the targets in the form p = k/2^n and tries to reduce the expression E(k/2^n) to an expression consisting of lower-level expressions:
E(1/32) = E(1/16) & R
E(2/32) = E(1/16) = E(2/16) & R
E(3/32) = E(3/16) & R
E(4/32) = E(2/16) = E(4/16) & R = E(2/2) & E(2/16)
E(5/32) = E(5/16) & R
E(6/32) = E(3/16) = E(6/16) & R = E(2/2) & E(3/16)
E(7/32) = E(7/16) & R
E(8/32) = E(4/16) = E(8/16) & R = E(2/2) & E(4/16)
E(9/32) = E(9/16) & R = E(3/4) & E(3/8)
E(10/32) = E(5/16) = E(10/16) & R = E(2/2) & E(5/16)
E(11/32) = E(11/16) & R
E(12/32) = E(6/16) = E(12/16) & R = E(2/2) & E(6/16)
E(13/32) = E(13/16) & R
E(14/32) = E(7/16) = E(14/16) & R = E(2/2) & E(7/16)
E(15/32) = E(15/16) & R = E(3/4) & E(5/8)
E(16/32) = E(8/16) = E(2/2) & E(8/16)
E(17/32) = E(1/4) | E(3/8)
E(18/32) = E(9/16) = E(2/2) & E(9/16)
E(19/32) = E(1/2) | E(3/16)
E(20/32) = E(10/16) = E(2/2) & E(10/16)
E(21/32) = E(3/4) & E(7/8) = E(1/2) | E(5/16)
E(22/32) = E(11/16) = E(2/2) & E(11/16)
E(23/32) = E(1/4) | E(5/8) = E(1/2) | E(7/16)
E(24/32) = E(12/16) = E(2/2) & E(12/16)
E(25/32) = E(1/4) | E(3/8) = E(1/2) | E(9/16)
E(26/32) = E(13/16) = E(2/2) & E(13/16)
E(27/32) = E(3/4) | E(3/8) = E(3/4) | E(3/8) = E(1/2) | E(11/16)
E(28/32) = E(14/16) = E(2/2) & E(14/16)
E(29/32) = E(1/4) | E(7/8) = E(1/2) | E(13/16) = E(3/4) | E(5/8)
E(30/32) = E(15/16) = E(2/2) & E(15/16)
E(31/32) = E(1/2) | E(15/16) = E(3/4) | E(7/8)
I hope it works, don't have more time to check it now. And I don't have a proof that it will always give results, but it seems to :-) |