# 
# 
# This script uses a simple recusrive algoritm to find G(N).  It uses
# a technique that looks at the factors of numbers in the sequence.
# We use the the basic observation that given a number N, if each prime 
# factor appears an even number of times, then N must be a perfect square.  
# For example, suppose N decomposes into prime factors f1, f2, ... fn, where 
# each factor appears an even number of times, designated by 2 * a1, 2 * a2, 
# etc:
# 
# N = (2 * a1) f1 * (2 * a2) f2 * ... * (2 * an) fn;
# 
# then:
# 
# N = [(a1 * f1)(a2 * f2)...(an * fn)] * [(a1 * f1)(a2 * f2)...(an * fn)];
# 
# making N a perfect square.
# 
# This algorithm works by factoring the initial N, noting which 
# factors appear an odd number of times, and adding numbers in the 
# sequence to "pair off" the odd factors.  Unfortunately, it is usally 
# the case that finding partners for some odd factors results in creating 
# new odd factors.  The algorthim does an exhaustive search of possble 
# partners to find the sequence with the smallest maximum value.
# 
# For example, consider G(5):
# 
# Non-trivial factors of 5 - (5)
# 
# The largest factor that appears an odd number of times is 5 - so we
# then add into the sequence the next highest number in sequence that
# provides 5 as a factor:
# 
# 5 * 6? (no 5 factor)
# 5 * 7? (no 5 factor)
# 5 * 8? (no 5 factor)
# 5 * 9? (no 5 factor)
# 5 * 10? (yes!)
# 
# Adding 10 effectively cancels out the 5, but now our factorization is
# 5 * 10 = 2 * 5 * 5 - so now we start over again looking for 2:
# 
# 5 * 10 * 6 (yes!)
# 
# Our factors now are: 2 * 2 * 3 * 5 * 5.  This takes care of our 2
# problem, but now there is a leftover 3.  So lets try again:
# 
# 5 * 6 * 10 * 6? (no, we already added 6)
# 5 * 6 * 10 * 7? (no 3 factor)
# 5 * 6 * 10 * 8? (no 3 factor)
# 5 * 6 * 10 * 9? (even number of 3, useless)
# 5 * 6 * 10 * 10? (no 3 factor)
# 5 * 6 * 10 * 11? (no 3 factor)
# 5 * 6 * 10 * 12? (yes)
# 
# Thus we now have 2 * 2 * 2 * 2 * 3 * 3 * 5 * 5 = (2 * 2 * 3 * 5) ^ 2 = 
# (60) ^2
# 
# We now note 12 as the "best" answer.  We will no longer consider
# values over 12.  Technically, the algorithm will check the next
# available integer for 3.  Because the next value is greater than our 
# "best" value of 12, we stop immediately:
# 
# 5 * 6 * 10 * 13? (stop, 13 > 12)
# 
# The algorithm will now backtrack to when it was checking "2"s:
# 
# 5 * 10 * 6 (backtrack through this)
# 5 * 10 * 7 (no 2)
# 5 * 10 * 8 (yes!)
# 
# Now we have 2 * 2 * 2 * 2 * 5 * 5 = (2 * 2 * 5) ^ 2 = (20) ^2
# 
# We now have an even better "best" G(N) of 10.  The algorithm will then 
# continue checking for 2's:
# 
# 5 * 10 * 9? (no 2)
# 5 * 10 * 10? (already had 10)
# 5 * 10 * 11 (no 11 > 10, end)
# 
# Now backtrack one last time, to when we just had 5:
# 
# 5 * 10? (backtrack through this)
# 5 * 11? (no, 11 > 10, so return)
# 
# Thus the algorithm ends with 10 as G(5).  Hopefully this explanation gives 
# the gist of what's going on.
# 
# --Ben Kennedy
# 

#
# here is the code 
#

use strict;

sub Graham {
  my $number = shift;

  # we will keep a hash of numbers that are in the sequence - this
  #  is so we can easily skip over numbers already in the sequence
  my %numbers = ( $number => 1 );

  # get list of factors in the number
  my @num_factors = @{get_factors($number)};

  # this hash will hold factors that appear an odd number of times
  my %factors;

  # first get the number of times each factor appears
  $factors{$_}++ foreach @{get_factors($number)}; 

  # now reduce it to simply a list of odd number facors
  %factors = map { ($_, 1) } grep { $factors{$_} % 2 } keys %factors;
  
  # record two things in state - the original starting value
  #  of the sequence, and a current "best" minimum G(N) - this
  #  will be populated after it finds the first sequence
  #
  my %state = ( 
	       CURRENT_MIN => '',
	       SEQ_START => $number
	      );
  
  # begin processing
  _rec_process(\%numbers, \%factors, \%state);

  # this sequence is not necessarily the shortest one
  # print join(' ', @{$state{CURRENT_MIN_SEQUENCE}}), "\n";

  # return the best minium value
  return $state{CURRENT_MIN};
}

sub _rec_process {
  my($numbers, $factors, $state) = @_;

  # optimization - only process each sequence of numbers once
  my $num_key = join(' ', sort { $a <=> $b } keys %$numbers);
  if ($state->{SEEN}{$num_key}++) {
    return;
  }

  # optimization - if sequence contains any values >= our
  #  best value, bail out.  This wont be processed if we
  #  haven't found a solution yet
  if ($state->{CURRENT_MIN}) {
    foreach my $num (keys %$numbers) {
      if ($num >= $state->{CURRENT_MIN}) {
	return;
      }
    }
  }

  # sort the factors by size - largest factor that needs a partner
  my @needed_facs = sort { $b <=> $a } keys %$factors;
  
  # if there are no needed factors, then we have a perfect square - note
  #  solution and return
  if (not @needed_facs) {
    my $new_best = (sort { $b <=> $a } keys %$numbers)[0];
    $state->{CURRENT_MIN} = $new_best;
    $state->{CURRENT_MIN_SEQUENCE} = [ sort { $a <=> $b } keys %$numbers ];
    return;
  } 

  # otherwise, try to find match for larger factors first.  Why larger?  
  #  This reduces run-time dramatically.  
  my $needed = shift(@needed_facs);

  # start at first multiple of needed factor after the start of the 
  #  sequence, and increase in multiples until we hit the highest 
  #  potential value.  If we don't have a current 'best' value, it 
  #  will not stop looking.  Fortunately, 4 * N will always be an upper 
  #  bound for this process.
  for (
       my $look_num = ($state->{SEQ_START} - 
		       ($state->{SEQ_START} % $needed) + $needed) ;
       not $state->{CURRENT_MIN} or ($look_num < $state->{CURRENT_MIN}) ;
       $look_num += $needed
      ) {
 
    # print "** $look_num - $needed\n";
   
    if ($numbers->{$look_num}) {
      # number is already in our sequence, move along
      next;
    }

    # now we count the number of times our needed number divides into
    #  $look_num - if it's an odd number of times, then it's a 
    #  suitable match
    #
    my $tmp = $look_num;
    my $count = 0;
    while(! ($tmp % $needed)) {
      $tmp /= $needed;
      $count++;
    }

    unless($count % 2) {
      # even number of times, so it won't help find a partner for the factor
      next;
    }
    
    # get all factors for this number so we can update the list
    my $look_factors = get_factors($look_num);

    # update factor list
    my %new_factors = %$factors;
    $new_factors{$_}++ foreach @$look_factors;      

    # remove all factors that are present an even number of times
    %new_factors = map { ($_, 1) } grep { $new_factors{$_} % 2 } 
      keys %new_factors;

    # recursively process - note new number and new list of odd
    #  factors
    _rec_process( { %$numbers, $look_num => 1 }, 
		  \%new_factors, $state);
  }
}


#
# Utility function that return the factors of a number as 
#  an array reference.  Caches the ref, so other code shouldn't 
#  alter it - this could be replaced with a CPAN module, but
#  factoring numbers of these sizes is relatively trivial
#

{
  my %factor_memo;
  sub get_factors {
    my $num = shift;    
    if ($factor_memo{$num}) {
      return $factor_memo{$num};
    }
    
    my @tmp;
    $factor_memo{$num} = \@tmp;
    rec_get_factors($num, \@tmp);
    return \@tmp;
  }
}

# yes, I like recursive algorithms
#
sub rec_get_factors {
  my($num, $ref) = @_;
  for (my $i = 2 ; $i <= $num ; $i++) {
    unless ($num % $i) {
      push @$ref, $i;
      rec_get_factors($num / $i, $ref);
      return;
    }
  }
}
