# Find a triangle-number decomposition by converting to # a Diophantine squares problem. use strict; # These have no quadratic residue for -1. # If they appear as factors in M an odd number of times # a solution is impossible. our @screeners = (3, 7, 11, 19, 23, 31); # The brute-force search for a quadratic solution. # We have reduced the number as much as we can before this # point. sub quad2 { my $M = shift; my $top = int(sqrt($M)); my $last = int($top/2) + 1; my ($i, $j, $rem); for ($i = $top; $i >= $last; --$i) { $rem = $M - $i*$i; $j = int(sqrt($rem)); if ($j*$j == $rem) { return ($i,$j); } } } # How many times does $n go into $M? sub countpowers { my ($M, $n) = @_; my $powercount = 0; while ($M % $n == 0) { $M = $M/$n; ++$powercount; } return ($M, $powercount); } # Reduce even numbers before solving by brute force. # Also, screen out some impossible cases. # Don't try a full factorization for now. sub quadreduce { my $M = shift; my $powercount; my $multfactor = 1; # Screen out odd impossible cases, and factor out pairs of 4k-1 factors. foreach my $num (@screeners) { ($M, $powercount) = countpowers($M, $num); if ($powercount % 2 == 1) { print "Eliminating impossible case with odd-power of 4k-1 factor $num\n"; return; } # In even cases, half the factors can be multiplied into the result. $powercount = $powercount/2; for (my $i = 0; $i < $powercount; $i++) { $multfactor *= $num; } } # Count powers of 2. We can handle the case of odd powers of 2. ($M, $powercount) = countpowers($M, 2); if ($powercount % 2 == 0) { # Two goes in an even number of times. my ($i, $j) = quad2($M); if ($powercount > 0) { $multfactor *= 1 << (($powercount)/2); } return ($i*$multfactor, $j*$multfactor); } # Two goes in an odd number of times. Use the i+j, i-j trick. my ($i, $j) = quad2($M); if ($powercount > 1) { $multfactor *= 1 << (($powercount - 1)/2); } return (($i + $j)*$multfactor, ($i - $j)*$multfactor); } my $M; $M = ($#ARGV >= 0)? $ARGV[0] : 987654321; # Convert from a triangle-number problem to a square-number problem. my $M2 = $M*8 + 3; # The top-level loop will try odd squares by brute force. my $top = int(sqrt($M2)); $top = $top - 1 if ($top % 2 == 0); # start with odd. my $last = int($top/2) + 1; my $k; for ($k = $top; $k >= $last; $k -= 2) { my $M3 = $M2 - $k*$k; my ($i, $j) = quadreduce($M3); if ($i) { # Convert solution back to triangle-numbers. my $new_i = ($i - 1)/2; my $new_j = ($j - 1)/2; my $new_k = ($k - 1)/2; print "Solution: triangle numbers $new_i, $new_j, $new_k\n"; my $tri_i = ($new_i+$new_i*$new_i)/2; my $tri_j = ($new_j+$new_j*$new_j)/2; my $tri_k = ($new_k+$new_k*$new_k)/2; my $sum = $tri_i + $tri_j + $tri_k; print "Verifying $tri_i + $tri_j + $tri_k = $sum = $M\n"; exit(0); } }