#!/usr/bin/perl
# The path where the binaries are:
$GGNFS_BIN_PATH="/afs/sipb.mit.edu/project/pari-gp/ggnfs/Linux/src";
# And some other popular choices:
#$GGNFS_BIN_PATH="../../src";
#$GGNFS_BIN_PATH="../ggnfs.vc/bin";
#$GGNFS_BIN_PATH="c:/mingw/msys/1.0/home/SamAdmin/ggnfs-0.73.1/src";
########################################################################
# factLat.pl
# Copyright 2004, Chris Monico.
#
#   This file is part of GGNFS.
#   GGNFS is free software; you can redistribute it and/or modify
#   it under the terms of the GNU General Public License as published by
#   the Free Software Foundation; either version 2 of the License, or
#   (at your option) any later version.
#
#   GGNFS is distributed in the hope that it will be useful,
#   but WITHOUT ANY WARRANTY; without even the implied warranty of
#   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
#   GNU General Public License for more details.
#
#   You should have received a copy of the GNU General Public License
#   along with GGNFS; if not, write to the Free Software
#   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
########################################################################
#  This script is known to work with Perl 5.8, and is known not to
#  work with 5.2. In between - I don't know.
########################################################################
use Math::BigInt;
use Math::BigFloat;
use Math::BigInt lib => 'GMP';

$SYS_BIN_PATH="c:/mingw/msys/1.0/bin";

$FORCECC="on"; # on | off | auto
$SAVEPAIRS=0;
$CLEANUP=0;
$PROMPTS=0;
$DOCLASSICAL=0;
$CHECK_BINARIES=1;
$ECHO_CMDLINE=1;

# If this is zero, the GGNFS polynomial selection code will be used (when needed).
# But the Kleinjung/Franke code is better, so you should use it if you can.
# i.e., try it - but if you seem to have some fatal errors, you can just change
# this to zero to revert to the GGNFS code.
$USE_KLEINJUNG_FRANKE_PS=1;

# Run the binaries at low priority?
#$NICE="nice -n 19 ";

# Set to 2 if you run into problems (which is definitely possible).
$LARGEP=3;

# The more extra FF's there are, the smaller the resulting pruned
# matrix will be (upto a certain point). Requiring the initial matrix
# to have an extra 10% columns is pretty reasonable.
$minextraFF=0.12;

# This should be in (0,1). Lower ==> smaller but denser pruned matrix
# (I think - this hasn't been extensively tested).
$matWtFactor=0.05;

# This controls how many partial relations can be combined to form full
# relations. The current built-in max of the software is 48, but this
# can result in very dense matrices. If you are getting matrices which
# are too dense (i.e., too dense to fit in RAM, or so dense that they
# are just taking a painfully long time to solve), you can drop this
# to, say 32, 24, or 20. The smaller this number, the more sieving you'll
# have to do, but you'll also get a sparser matrix.
$maxRelsInFF=28;

# This is for an Athlon 2800+ laptop. If your machine is about half as fast,
# replace this with a 2. 25% as fast, replace with a 4. It controls how long
# the polynomial selection phase will last.
$polySelTimeMultiplier=1.0;

################################################################
# Nothing configurable below here - don't mess with it unless  #
# you're fixing a bug or adding functionality.                 #
################################################################

if($^O ne "MSWin32") {
  $CAT="cat";
  $GZIP="gzip";
  $EXEC_SUFFIX="";
}
else {
  $CAT=$SYS_BIN_PATH."/cat.exe";
  $GZIP=$SYS_BIN_PATH."/gzip.exe";
  $EXEC_SUFFIX=".exe";
}

$LATSIEVER_L1=$GGNFS_BIN_PATH."/gnfs-lasieve4I12e".$EXEC_SUFFIX;
$LATSIEVER_L2=$GGNFS_BIN_PATH."/gnfs-lasieve4I13e".$EXEC_SUFFIX;
$LATSIEVER_L3=$GGNFS_BIN_PATH."/gnfs-lasieve4I14e".$EXEC_SUFFIX;
$LATSIEVER=$LATSIEVER_L1; # Just a default.

####################################################################
$MAKEFB=$GGNFS_BIN_PATH."/makefb".$EXEC_SUFFIX;
$PROCRELS=$GGNFS_BIN_PATH."/procrels".$EXEC_SUFFIX;
$CLSIEVE=$GGNFS_BIN_PATH."/sieve".$EXEC_SUFFIX;
$MATBUILD=$GGNFS_BIN_PATH."/matbuild".$EXEC_SUFFIX;
$MATSOLVE=$GGNFS_BIN_PATH."/matsolve".$EXEC_SUFFIX;
$SQRT=$GGNFS_BIN_PATH."/sqrt".$EXEC_SUFFIX;
$POLYSELECT=$GGNFS_BIN_PATH."/polyselect".$EXEC_SUFFIX;
$POL51M0=$GGNFS_BIN_PATH."/pol51m0b".$EXEC_SUFFIX;
$POL51OPT=$GGNFS_BIN_PATH."/pol51opt".$EXEC_SUFFIX;
$PLOT=$GGNFS_BIN_PATH."/autogplot.sh";
$DEFAULT_PAR_FILE=$GGNFS_BIN_PATH."/def-par.txt";
$DEFAULT_POLSEL_PAR_FILE=$GGNFS_BIN_PATH."/def-nm-params.txt";

$DEPFILE="deps";
$SPMAT="spmat";
$RELSBIN="rels.bin";
$LOGFILE="ggnfs.log";
$LARGEPRIMES="-".$LARGEP."p";

# This file is used by this script to detect
# changes in parameter settings.
$PARAMFILE=".params";

# This is to tell the lattice siever where to dump the next special-q
# if it's interrupted, for example, with CTRL-C.
$PNUM=0;


                                                                                                
###########################################################
sub gcd {
###########################################################
# This will be used to make a function which peels off prime factors
# when the square root step returns composite factors, so that
# (1) We can report the prime factors themselves.
# (2) We know when enough square roots have been computed.
###########################################################
# Sample usage: (Note that the declarations must be this way,
# to force the variables to have the right type!).
#$x= Math::BigInt->new($ARGV[0]);
#$y= Math::BigInt->new($ARGV[1]);
#$g=gcd($x, $y);
###########################################################
  if ($_[0]==0) { return $_[1]; }
  if ($_[1]==0) { return $_[0]; }
  return gcd($_[1]%$_[0], $_[0]);
}

#-------------------------------------------------------------------------------
#probab_prime_p(n, reps)
#  Return whether n is probably prime or not.
#    '' = composite
#    1 = probably prime
#    2 = prime
#  This function uses GMP PERL MODULE (gmp-4.1.4/demos/perl) if available.
#  GMP::Mpz::probab_prime_p is very fast but never returns 2 because it's a bool
#  function.
if (eval('use GMP::Mpz;1;')) {
  *probab_prime_p = sub {
    GMP::Mpz::probab_prime_p("$_[0]", $_[1]);  #convert Math::BigInt to scalar
  };
} else {
  *probab_prime_p = sub {
    my ($n, $reps) = @_;
    ref $n or $n = new Math::BigInt($n);
    $n->is_negative() and $n->babs();
    #Trial division
    if ($n < 1000000) {
      $n = "$n" + 0;  #scalar
      $n < 2 and return '';
      foreach (@primes1000) {
        my $q = int($n / $_);
        $q < $_ and return 2;
        $n - $_ * $q or return '';
      }
      return 2;
    }
    $n->bgcd($primorial168)->is_one() or return '';  #bgcd returns new object
    #Fermat test
    #  gcd(n,210)==1 is ensured by trial division.
    my $nm1 = $n->copy()->bdec();  #new object
    new Math::BigInt(210)->bmodpow($nm1, $n)->is_one() or return '';
    #Miller Rabin test
    #  Reference: gmp-4.1.4/mpz/millerrabin.c
    my $k = scan1($nm1, 0);
    my $q = $nm1->copy()->brsft($k, 2);  #new object
  MILLER_RABIN_LOOP:
    while ($reps-- > 0) {
      my $x;
      do {
        $x = urandomb(sizeinbase($n, 2) - 1);
      } while ($x <= 1);
      $x->bmodpow($q, $n)->is_one() || $x == $nm1 and next MILLER_RABIN_LOOP;
      my $i;
      for ($i = 1; $i < $k; $i++) {
        $x->bmul($x)->bmod($n) == $nm1 and next MILLER_RABIN_LOOP;
        $x->is_one() and return '';
      }
      return '';
    }
    1;
  };
  #primes1000
  #  Primes up to 1000.
  @primes1000 = (
    2, 3, 5, 7, 11, 13, 17, 19, 23, 29,
    31, 37, 41, 43, 47, 53, 59, 61, 67, 71,
    73, 79, 83, 89, 97, 101, 103, 107, 109, 113,
    127, 131, 137, 139, 149, 151, 157, 163, 167, 173,
    179, 181, 191, 193, 197, 199, 211, 223, 227, 229,
    233, 239, 241, 251, 257, 263, 269, 271, 277, 281,
    283, 293, 307, 311, 313, 317, 331, 337, 347, 349,
    353, 359, 367, 373, 379, 383, 389, 397, 401, 409,
    419, 421, 431, 433, 439, 443, 449, 457, 461, 463,
    467, 479, 487, 491, 499, 503, 509, 521, 523, 541,
    547, 557, 563, 569, 571, 577, 587, 593, 599, 601,
    607, 613, 617, 619, 631, 641, 643, 647, 653, 659,
    661, 673, 677, 683, 691, 701, 709, 719, 727, 733,
    739, 743, 751, 757, 761, 769, 773, 787, 797, 809,
    811, 821, 823, 827, 829, 839, 853, 857, 859, 863,
    877, 881, 883, 887, 907, 911, 919, 929, 937, 941,
    947, 953, 967, 971, 977, 983, 991, 997,
  );
  #primorial168
  #  Product of primes up to 1000.
  $primorial168 = Math::BigInt->bone();
  foreach (@primes1000) {
    $primorial168->bmul($_);
  }
  #scan1(n, start)
  #  Return the index of the least significant 1 in base 2.
  *scan1 = sub {
    my ($n, $start) = @_;
    $n = new Math::BigInt($n);  #new object
    $n->is_negative() and $n->babs();
    $start > 0 and $n->brsft($start, 2);
    $n->is_zero() and return 0xffffffff;
    my $q;
    ($n, $q) = $n->bdiv(0x100000000);
    until ($q) {
      $start += 32;
      ($n, $q) = $n->bdiv(0x100000000);
    }
    until ($q & 1) {
      $start++;
      $q >>= 1;
    }
    $start;
  };
  #sizeinbase(n, base)
  #  Return the size of n measured in number of digits in base base.
  *sizeinbase = sub {
    my ($n, $base) = @_;
    $n = new Math::BigInt($n);  #new object
    $n->is_negative() and $n->babs();
    $base < 2 and $base = 2;
    $base = new Math::BigInt($base);  #new object
    my @list = ();
    while ($base <= $n) {
      push(@list, $base);
      $base *= $base;  #new object
    }
    my $size = 1;
    while (@list) {
      $base = pop(@list);
      if ($base <= $n) {
        $n->bdiv($base);
        $size += 1 << @list;
      }
    }
    $size;
  };
  #urandomb(size)
  #  Generate a random integer in the range 0 to 2^size-1, inclusive.
  *urandomb = sub {
    my ($size) = @_;
    my $n = new Math::BigInt(int(rand(1 << ($size & 31))));
    while ($size >= 32) {
      $n->bmul(0x100000000)->bior(int(rand(0x100000000)));  #don't use blsft
      $size -= 32;
    }
    $n;
  };
}
#-------------------------------------------------------------------------------
###########################################################
sub getPrimes {
###########################################################
# Read the log file to see if we've found all the prime
# divisors of N yet.
###########################################################

  open(INFO,$LOGFILE);
  while (<INFO>) {
    chomp;
    @_ = split;
    if (/r\d=/) {
      s/.*r\d=//;
      if ((length($_) > 1) && (length($_) < length($N))) {
        # Is this a prime divisor or composite?
        if (/pp/) {
          s/\(.*\)//; # Strip off the (pp <digits>) part.
          s/\s*//g; # Remove whitespace.
          # If this is a prime we don't already have, add it.
          my $found=0;
          for ($i=0; $i<=$#PRIMES; $i++) {
            if ($_ == $PRIMES[$i]) { $found=1; }
          }
          if (!($found)) { push(@PRIMES, $_); }
        } else {
          s/\(.*\)//; # Strip off the (c<digits>) part.
          s/\s*//g; # Remove whitespace.
          push(@COMPS, $_);
        }
      }
    }
  }
  close(INFO);
  # Now, try to figure out if we have all the prime factors:
  my $x= Math::BigInt->new('1');
  for ($i=0; $i<=$#PRIMES; $i++) {
    $x = $x*$PRIMES[$i];
  }
  if ($x==$N || probab_prime_p($N/$x, 10)) { 
    $x==$N or push(@PRIMES, $N/$x);
    open(OF, ">>$LOGFILE");
    while ($_ = shift @PRIMES) {
      printf(OF "-> p: $_ (pp%d)\n", length($_));
    }
    return 1; 
  }
  # Here, we could try to recover other factors by division,
  # but until we have a primality test available, this would
  # be pointless since we couldn't really know if we're done.
  return 0; 
    
}

###########################################################
sub sigDie {
###########################################################
  die "Signal caught. Terminating...\n";
}

######################################################
sub loadDefaultParams {
######################################################
# 
# These are default parameters for different size factorizations.
# They will be used only to fill in any non-user-supplied parameters.
# The format of the file is as follows:
# type, digits, deg, maxs1, maxskew, goodScore, eFrac, j0, j1, eStepSize, maxTime,
#       rlim, alim, lpbr, lpba, mfbr, mfba, rlambda, alambda, qintsize
# where 'type' is gnfs or snfs.
# arg0 = number of digits in N.
# arg1 = type

  die("loadDefaultParams(): Insufficient arguments!\n") 
      unless ($#_ >=1);
  my $DIGS=$_[0];
  my $type=$_[1];
  die("Could not find default parameter file $DEFAULT_PAR_FILE!\n") 
      unless (-e $DEFAULT_PAR_FILE);
  open(IF, $DEFAULT_PAR_FILE);
  my $howClose=1000;
  while (<IF>) {
    s/#.*//; # Remove comments
    s/\s*//g; # Remove whitespace.
    if (length($_)>0) {
      @_ = split /,/;
      my $t=$_[0];
      if ($t eq $type) {
        my $d=$_[1];
        if (abs($d-$DIGS)<$howClose) {
          $o=2;
          $howClose=abs($d-$DIGS);
          $digLevel=$d;
          $DEG=$_[$o+0]; 
          $MAXS1=$_[$o+1]; 
          $MAXSKEW=$_[$o+2];
          $GOODSCORE=$_[$o+3];  
          $EFRAC=$_[$o+4];
          $J0=$_[$o+5]; $J1=$_[$o+6];
          $ESTEPSIZE=$_[$o+7]; 
          $MAXTIME=$_[$o+8];
          $RLIM=$_[$o+9];     $ALIM=$_[$o+10]; 
          $LPBR=$_[$o+11];    $LPBA=$_[$o+12];
          $MFBR=$_[$o+13];    $MFBA=$_[$o+14];
          $RLAMBDA=$_[$o+15]; $ALAMBDA=$_[$o+16]; 
          $QINTSIZE=$_[$o+17];
          $classicalA=$_[$o+18];
          $classicalB=$_[$o+19];
          $QSTEP=$QINTSIZE;
        }
      }
    }
  }
  close(IF);
  printf "-> Selected default factorization parameters for $digLevel digit level.\n";
  if ($type eq "gnfs") {
    if ($DIGS < 110) { $LATSIEVER=$LATSIEVER_L1; }
    elsif ($DIGS < 135) { $LATSIEVER=$LATSIEVER_L2; }
    else { $LATSIEVER=$LATSIEVER_L3; }
  } else {
    if ($DIGS < 150) { $LATSIEVER=$LATSIEVER_L1; }
    elsif ($DIGS < 180) { $LATSIEVER=$LATSIEVER_L2; }
    else { $LATSIEVER=$LATSIEVER_L3; }
  }
  printf "-> Selected lattice siever: $LATSIEVER\n";
}

######################################################
sub loadPolselParamsPol5 {
######################################################
# 
# These are default parameters for polynomial selection using the
# Kleinjung/Franke tool. 
# arg0 = number of digits in N.

  die("loadDefaultParams(): Insufficient arguments!\n") 
      unless ($#_ >=0);
  my $DIGS=$_[0];
  die("Could not find default parameter file $DEFAULT_POLSEL_PAR_FILE!\n") 
      unless (-e $DEFAULT_POLSEL_PAR_FILE);
  open(IF, $DEFAULT_POLSEL_PAR_FILE);
  my $howClose=1000;
  while (<IF>) {
    s/#.*//; # Remove comments
    s/\s*//g; # Remove whitespace.
    if (length($_)>0) {
      @_ = split /,/;
      my $d=$_[0];
      if (abs($d-$DIGS)<$howClose) {
        $o=1;
        $howClose=abs($d-$DIGS);
        $digLevel=$d;
        $maxPSTime=60*$_[$o+0];
        $search_a5step=$_[$o+1]; 
        $npr=$_[$o+2]; 
        $normmax=$_[$o+3]; 
        $normmax1=$_[$o+4];
        $normmax2=$_[$o+5];
        $murphymax=$_[$o+6];
      }
    }
  }
  close(IF);
  printf "-> Selected default polsel parameters for $digLevel digit level.\n";
}

#######################################################
sub terminate_search {
    $terminate_job=1;
    printf "Terminated on $pname by SIGTERM\n";
}

#######################################################
sub runPol5 {
  $projectname=$NAME.".polsel";
  use Sys::Hostname;
  $host=hostname;
  $pname="$projectname.$host.$$";

  open(OF, ">$pname.data");
  printf OF "N ".$N;
  close(OF);
  my $terminate_job=0;
  local $SIG{'TERM'}='terminate_search';

  loadPolselParamsPol5(length($N));
  $maxPSTime *= $polySelTimeMultiplier;
  my $hmult=1e3;
  loadDefaultParams(length($N), "gnfs");
  my %bestpolyinf = ();
  $bestpolyinf{Murphy_E} = 0;

  my $H=0;
  my $startTime = time;
  for(my $nerr=0;$terminate_job==0 && $nerr<2;) {
    my $HH=$H+$search_a5step;
    printf "-> Searching leading coefficients from %d to %d.\n", $H*$hmult+1, $HH*$hmult;
    $cmd="$NICE \"$POL51M0\" -b $pname -v -v -p $npr -n $normmax -a $H -A $HH > $pname.log";
    printf("=> $cmd\n");
    my $res=system($cmd);
    die "Return value $res. Terminating...\n" if ($res);
    $nerr=0;
    open(GR,"$pname.log");
    my @logout=<GR>;
    close(GR);
    $suc = grep(/success/, @logout);
    if($suc!=0) {
      $cmd="$NICE \"$POL51OPT\" -b $pname -v -v -n $normmax1 -N $normmax2 -e $murphymax > $pname.log";
      printf("=> $cmd\n");
      $res=system($cmd);
      die "Return value $res. Terminating...\n" if ($res);
      system("\"$CAT\" $pname.51.m >> $projectname.51.m.all");
      open(GR,"<$pname.cand");
      my %polyinf = ();
      my $changed = 0;
      while(<GR>) {
        chomp; s/\r?$//;
        if (s/^BEGIN POLY//) {
          s/^ #//;
          %polyinf = split;
        } elsif (/^END POLY/) {
          if ($polyinf{Murphy_E} > $bestpolyinf{Murphy_E}) {
            %bestpolyinf = %polyinf;
            $changed = 1;
          }
          %polyinf = ();
        } else {
          my ($key, $val) = split;
          $key =~ s/^X/c/;
          $polyinf{$key} = $val;
        }
      }
      close(GR);
      if ($changed) {
        foreach my $key (sort keys %bestpolyinf) {
          print "$key: $bestpolyinf{$key}\n";
        }
        open(BP, ">$NAME.poly");
        print BP "name: $NAME\n";
        print BP "n: $N\n";
        foreach my $key (reverse sort keys %bestpolyinf) {
          if ($key =~ /c\d+/ or $key =~ /Y\d+/) {
            print BP "$key: $bestpolyinf{$key}\n";
          } elsif ($key =~ /skewness/) {
            print BP "skew: $bestpolyinf{$key}\n";
          } else {
            print BP "# $key $bestpolyinf{$key}\n";
          }
        } 
        print BP "type: gnfs\n";
        print BP "rlim: $RLIM\n";
        print BP "alim: $ALIM\n";
        print BP "lpbr: $LPBR\n";
        print BP "lpba: $LPBA\n";
        print BP "mfbr: $MFBR\n";
        print BP "mfba: $MFBA\n";
        print BP "rlambda: $RLAMBDA\n";
        print BP "alambda: $ALAMBDA\n";
        print BP "qintsize: $QINTSIZE\n";
        close(BP);
      }
      system "\"$CAT\" $pname.cand >> $projectname.cand.all";
      unlink "$pname.cand";
    }
    printf "-> =====================================================\n";
    printf("-> Best score so far: %e (goodScore=%e)\n",$bestpolyinf{Murphy_E},$murphymax);
    printf "-> =====================================================\n\n";
    unlink "$pname.log";
    unlink "$pname.51.m";
    my $nowTime = time;
    if ($nowTime - $startTime > $maxPSTime) {
      $terminate_job=1;
    }
    $H=$HH;
  }
  unlink "$pname.data";

  # What remains to be done is to:
  # (1) Search the x.cand.x file for the best candidate.
  # (2) load default factorization parameters.
  # (3) Output a GGNFS polynomial file with the poly and parameters.
  # (4) Delete the intermediate files.
  # In fact, we should probably do (1) above and some file renaming,
  # so that we don't have to search through the whole (growing) file
  # at each iteration. Keep the best poly in a seperate file and
  # maybe even just delete the other candidates.
  #
  # S.Chong: 
  # (1) through (3) are done, (4) if you don't care about leaving the *.all
  # files lying around.  Only the best candidate is saved in a .poly file,
  # to look at the rest you'll have to dig through the *.cand.all file.
  #     The next thing to do is figure out how to support multiple polysel
  # clients.  The filenames shouldn't need to be changed, but we'll need to
  # get/update $H from a lock-protected file and also lock-protect the *.all
  # files.  Then the rest should be easy.  I imagine multi-client sieving
  # could be enhanced similarly around $Q0.
}

######################################################
sub runPolyselect {
######################################################
# We will start with a higher leading coefficient divisor. When it
# appears that we are searching in an interesting range, it will
# be backed down so that the resulting range can be searched with
# a finer resolution. This means that from time to time, the same
# poly will be found several times as we hone in on a region.
  my @lcdChoices=(2,4,4,12,12,24,24,48,48,144,144,720,5040);
  my $lcdLevel=4+(length($N)-70)/10;
  if ($lcdLevel < 0) { $lcdLevel=0;}
  if ($lcdLevel > 12) { $lcdLevel=12;}
  my $E0=1;
  my $firstGoodTime=0;
  printf "-> Starting search with leading coefficient divisor $lcdChoices[$lcdLevel].\n";

  loadDefaultParams(length($N), "gnfs");
  $MAXTIME *= $polySelTimeMultiplier;
  my $E1=$E0+$ESTEPSIZE;
  my $goodPFound=0;
  my $bestScore=0.0;
  my $startTime = time;
  my $done=0;
  my $bestLC=1;
  my $multiplier=0.75;
  while (!$done) {
    $LCD=$lcdChoices[$lcdLevel];
    open(OF, ">$NAME.polsel");
    printf OF "name: $NAME\n";
    printf OF "n: $N\n";
    printf OF "deg: $DEG\n";
    printf OF "bf: best.poly\n";
    printf OF "maxs1: $MAXS1\n";
    printf OF "maxskew: $MAXSKEW\n";
    printf OF "enum: $LCD\n";
    printf OF "e0: $E0\n";
    printf OF "e1: $E1\n";
    $CUTOFF=0.75*$GOODSCORE;
    printf OF "cutoff: $CUTOFF\n";
    printf OF "examinefrac: $EFRAC\n";
    printf OF "j0: $J0\n";
    printf OF "j1: $J1\n";
    close(OF);
    my $cmd="$NICE \"$POLYSELECT\" -if $NAME.polsel";
    print "=>$cmd\n" if($ECHO_CMDLINE);
    $res=system($cmd);
    die "Return value $res. Terminating...\n" if ($res);
    # Find the score of the best polynomial:
    # E(F1,F2) = 
    open(INFO, "best.poly");
    my @polyinf=<INFO>;
    close(INFO);
    my @TMP=grep(/E\(F1,F2\) =/, @polyinf);
    my $SCORE=$TMP[0];
    $SCORE =~ s/.*E\(F1,F2\) =//;
    @TMP=grep(/^c$DEG/, @polyinf);
    $_=$TMP[0];
    /^c$DEG:\s(\d*)/;
    $LC = $1;
    if (($SCORE > $multiplier*$GOODSCORE)&&($LC != $bestLC) && ($LC != $lastLC)) {
      $multipler *= 1.1;
      if ($multiplier > 0.9) { $multiplier = 0.9; }
      $lastLC = $LC;
      my $newLCDLevel = $lcdLevel-1;
      if ($newLCDLevel < 0) { $newLCDLevel=0; }
      $E0 = $E0*$lcdChoices[$lcdLevel]/$lcdChoices[$newLCDLevel] - $ESTEPSIZE;
      if ($lcdLevel != $newLCDLevel) {
        printf "-> Leading coefficient divisor dropped from $lcdChoices[$lcdLevel] to $lcdChoices[$newLCDLevel].\n";
      }
      $lcdLevel = $newLCDLevel;
    }
    if (($SCORE > $bestScore)&&($LC != $bestLC)) { 
      $bestScore = $SCORE;
      $bestLC=$LC;
      $lastBestTime = time;
      rename "best.poly", "thebest.poly";
      # We should now fill in the missing parameters with the defaults
      # loaded in from table. Do this now, so that the user has the
      # option to kill the script and still have a viable poly file.
      open(IF, "thebest.poly");
      open(OF, ">$NAME.poly");
      while (<IF>) {
        chomp;
        if (/rlim:/) { $_="rlim: $RLIM"; }
        if (/alim:/) { $_="alim: $ALIM"; }
        if (/lpbr:/) { $_="lpbr: $LPBR"; }
        if (/lpba:/) { $_="lpba: $LPBA"; }
        if (/mfbr:/) { $_="mfbr: $MFBR"; }
        if (/mfba:/) { $_="mfba: $MFBA"; }
        if (/rlambda:/) { $_="rlambda: $RLAMBDA"; }
        if (/alambda:/) { $_="alambda: $ALAMBDA"; }
        if (/qintsize:/) { $_="qintsize: $QINTSIZE"; }
        printf OF "$_\n";
      }
      printf OF "type: gnfs\n";
      close(OF);
      close(IF);
      unlink "thebest.poly";
      unlink "best.poly";
    }
    printf "-> =====================================================\n";
    printf("-> Best score so far: %f (goodScore=%f)     \n",$bestScore,$GOODSCORE);
    printf "-> =====================================================\n";

    $E0 += $ESTEPSIZE;
    $E1 = $E0 + $ESTEPSIZE;
    $done=0;
    if ($bestScore > 1.4*$GOODSCORE) { $goodPFound=1; }
    if ($goodPFound) {
      # We will allow another 5 minutes just in case there happens to
      # be a really good poly nearby (or the 'goodScore' value was too low)
      my $elapsed = time - $lastBestTime;
      if ($elapsed > 300) { 
        $done=1;
      }
    }
    my $nowTime=time;
    if ($nowTime - $startTime > $MAXTIME) { $done=1; }
  }
  printf("-> Using poly with score=%f\n", $bestScore);

  unlink "$NAME.polsel";
}

######################################################
sub changeParams {
######################################################
# Change the factor base sizes for a factorization
# which has already been started. This is done by
# dumping the raw siever output from the rels.bin files,
# deleting the rels.bin files, re-setting the parameters,
# and reprocessing all of the relations.
  open(LF, ">>$LOGFILE");
  print LF "Parameters change detected.\n";
  close(LF);

  print "-> Parameter change detected...\n";
  print "-> Dumping relations...\n";
  $cmd="$NICE \"$PROCRELS\" -fb $NAME.fb -prel $RELSBIN -dump";
  print "=>$cmd\n" if($ECHO_CMDLINE);
  $res=system($cmd);
  die "Return value $res. Terminating...\n" if ($res);
  unlink <$RELSBIN*>, <cols*>, <deps*>, 'factor.easy', <lpindex*>;
  unlink <$NAME.*.afb.*>;

  print "-> Making new factor base files...\n";
  $cmd="$NICE \"$MAKEFB\" -rl $RLIM -al $ALIM -lpbr $LPBR -lpba $LPBA $LARGEPRIMES -of $NAME.fb -if $NAME.poly";
  print "=>$cmd\n" if($ECHO_CMDLINE);
  $res=system($cmd);
  die "Return value $res. Terminating...\n" if ($res);

  print "-> Reprocessing siever output...\n";
  $i=0;
  while (-e "spairs.dump.$i") {
    $cmd="$NICE \"$PROCRELS\" -fb $NAME.fb -prel $RELSBIN -newrel spairs.dump.$i";
    print "=>$cmd\n" if($ECHO_CMDLINE);
    $res=system($cmd);
    $i++;
  }
  # Update the paramfile, used by this script to detect changes
  # in parameter settings.
  open(INFO, ">$PARAMFILE");
  print INFO "rlim: $RLIM\n"; 
  print INFO "alim: $ALIM\n"; 
  print INFO "lpbr: $LPBR\n"; 
  print INFO "lpba: $LPBA\n"; 
  close(INFO);
} 
  
###########################################################
sub checkParamFile {
###########################################################
# Check two things:
# (1) Does a $PARAMFILE exist? If not, create it.
# (2) Are the values in $PARAMFILE the same as the current
#     values? If not, call changeParams to sync everything up.
  if (!(-e $PARAMFILE)) {
    print "-> Creating param file to detect parameter changes...\n";
    open(INFO, ">$PARAMFILE");
    print INFO "rlim: $RLIM\n"; 
    print INFO "alim: $ALIM\n"; 
    print INFO "lpbr: $LPBR\n"; 
    print INFO "lpba: $LPBA\n"; 
    close(INFO);
    return ;
  }
  # Okay - it exists. We should check it to see if the
  # parameters are the same as the current ones.
  open(INFO, $PARAMFILE);
  my @lastParams=<INFO>;
  close(INFO);
  my @TMP=grep(/rlim:/, @lastParams);
  my $OLDRLIM=$TMP[0];
  $OLDRLIM =~ s/.*rlim://;
  @TMP=grep(/alim:/, @lastParams);
  my $OLDALIM=$TMP[0];
  $OLDALIM =~ s/.*alim://;
  @TMP=grep(/lpbr:/, @lastParams);
  my $OLDLPBR=$TMP[0];
  $OLDLPBR =~ s/.*lpbr://;
  @TMP=grep(/lpba:/, @lastParams);
  my $OLDLPBA=$TMP[0];
  $OLDLPBA =~ s/.*lpba://;
  if ($OLDRLIM != $RLIM || $OLDALIM != $ALIM ||
      $OLDLPBR != $LPBR || $OLDLPBA != $LPBA) {
    changeParams;
  } else {
    print "-> No parameter change detected. Resuming.\n";
  }
}

###########################################################
sub plotLP {
###########################################################
  if ($CHECK_BINARIES) {
    return unless (-x $PLOT);
  }
  open(OUTF, ">.lprels") || return;
  print OUTF "0, 0\n";
  open(OUTRELS, ">.rels") || return;
  print OUTRELS "0, 0\n";
  open(LOG, $LOGFILE) || return;
  while ($_ = <LOG>) {
    chomp;
    tr/a-z/A-Z/; # Convert to upper case.
    if ( /LARGEPRIMES/) {
      s/\[.*]\s*//;  # Strip out the date.
      s/(LARGEPRIMES: |RELATIONS: )//g;  # Strip out the labels.
      s/\s//g; # Remove whitespace.
      @_ = split /,/;
      $excessLP=$_[0]-$_[1];
      print OUTF "$_[1], $excessLP\n";
    } elsif (/FINALFF/) {
      s/\[.*]\s*//;  # Strip out the date.
      s/(RELS:|INITIALFF:\d*,|FINALFF:)//g;  # Strip out the labels.
      s/\s//g; # Remove whitespace.
      @_ = split /,/;
      print OUTRELS "$_[0], $_[1]\n";
    }
  }
  close(OUTF);
  close(OUTRELS);
  close(LOG);
  $ENV{'XAXIS'}='Total relations';
  $ENV{'YAXIS'}="";
  $cmd="\"$PLOT\""." xprimes.jpg 'ExcessLargePrimes' .lprels";
  print "=>$cmd\n" if($ECHO_CMDLINE);
  $res=system($cmd);
  die "Return value $res. Terminating...\n" if ($res);
  $ENV{'XAXIS'}='Total relations';
  $ENV{'YAXIS'}="Full relation-sets";
  $cmd="\"$PLOT\""." relations.jpg 'TotalFF' .rels";
  print "=>$cmd\n" if($ECHO_CMDLINE);
  $res=system($cmd);
  die "Return value $res. Terminating...\n" if ($res);
  unlink '.lprels';
  unlink '.rels';
}

######################################################
sub checkParams {
######################################################
  if ($CHECK_BINARIES) {
    $missing .= 'makefb ' unless (-x $MAKEFB);
    $missing .= 'procrels ' unless (-x $PROCRELS);
    $missing .= 'matbuild ' unless (-x $MATBUILD);
    $missing .= 'matsolve ' unless (-x $MATSOLVE);
    $missing .= 'sqrt ' unless (-x $SQRT);
    $missing .= '(lattice siever) ' unless (-x $LATSIEVER);
  }
  if ($missing) {
    print "-> Could not find GGNFS programs: $missing.\n";
    print "-> Did you set GGNFS_BIN_PATH properly in this script?\n";
    print "-> It is currently set to:";
    print "-> GGNFS_BIN_PATH=$GGNFS_BIN_PATH\n";
    exit -1;
  }
  die("Error: 'n' not supplied!\n") unless ($N);
  die("Error: 'm' not supplied!\n") unless ($M || defined $COEFHASH{Y1});
  die("Error: polynomial not supplied!\n") unless (defined $COEFHASH{'c'.$DEGREE});
  die("Error: 'skew' not supplied!\n") unless ($SKEW);
  die("Error: 'rlim' not supplied!\n") unless ($RLIM);
  die("Error: 'alim' not supplied!\n") unless ($ALIM);
  die("Error: 'lpbr' not supplied!\n") unless ($LPBR);
  die("Error: 'lpba' not supplied!\n") unless ($LPBA);
  die("Error: 'mfbr' not supplied!\n") unless ($MFBR);
  die("Error: 'mfba' not supplied!\n") unless ($MFBA);
  die("Error: 'rlambda' not supplied!\n") unless ($RLAMBDA);
  die("Error: 'alambda' not supplied!\n") unless ($ALAMBDA);
  die("Error: 'qintsize' not supplied!\n") unless ($QSTEP);
}

######################################################
sub makeJobFile {
######################################################
# arg0 = filename
# arg1 = first q value
# arg2 = q range size
# arg3 = client num
# arg4 = number of clients
# arg5 = A0 value for classical sieving
# arg6 = A1 value for classical sieving
# arg7 = B0 value for classical sieving
# arg8 = B1 value for classical sieving
# Note: arguments 5-8 are only used if arg1=arg2=0. Otherwise,
# they need not even be supplied.
  if ($#_ < 4) { die("makeJobFile() : Not enough arguments!\n") };
  $FNAME=$_[0];
  my $firstQ = $_[1];
  my $qrSize = $_[2];
  my $client = $_[3]-1; # Convert to 0,1,..., numClients-1.
  my $numClients = $_[4];
  my $A0=0; my $A1=0;
  my $B0=0; my $B1=0;
  if ($#_ >= 8) {
    $A0 = $_[5]; $A1 = $_[6];
    $B0 = $_[7]; $B1 = $_[8];
  }
  my $sieveType=0;  

  if (($firstQ >0) || ($qrSize > 0)) {
    $sieveType=1;
    # First, find the proper q0 and qrSize for this client.
    # The idea is that client 'c' should sieve over ranges
    #     [QSTART + k*qrSize, QSTART + (k+1)qrSize]
    # with k == c (mod numClients). Thus, we need to find the
    # first such range containing q0.
    my $k=int(($firstQ - $QSTART)/$qrSize);
    if ($k*$qrSize + $QSTART > $firstQ) {
      # Could this happen from rounding error? 
      $k--;
    }
    while (($k % $numClients) != $client) {
      $k += 1;
    }
    $q0 = $QSTART + $k*$qrSize;
    $q1 = $q0 + $qrSize;
    printf "-> makeJobFile(): q0=$q0, q1=$q1.\n";
    if ($firstQ >= $q1) {
      $k += $numClients;
      $q0 = $QSTART + $k*$qrSize;
      $q1 = $q0 + $qrSize;
    } else {
      if ($firstQ > $q0) {
        $q0 = $firstQ;
      }
      $qrSize = $q1 -$q0;
    }
    printf "-> makeJobFile(): Adjusted to q0=$q0, q1=$q1.\n";
    $Q0 = $q0;
    $Q1 = $q1;
    $thisQRSize = $qrSize;
  }

  unlink $FNAME;
  open(OUTF, ">$FNAME");
  print OUTF "n: $N\nm: $M\n";
  # The polynomial coefficients:
  for ($i=0; length(${COEF[$i]}) > 0; $i += 2) {
    print OUTF "${COEF[$i]} ${COEF[${i}+1]}\n";
  }
  print OUTF "skew: $SKEW\n";
  print OUTF "rlim: $RLIM\n";

  if (($sieveType==1) && ($ALIM > $q0)) { 
    $sieverAL = $q0-1; 
    unlink <$JOBNAME.afb.*>;
  }
  else { $sieverAL = $ALIM; }
  print OUTF "alim: $sieverAL\n";
  print OUTF "lpbr: $LPBR\n";
  print OUTF "lpba: $LPBA\n";
  print OUTF "mfbr: $MFBR\n";
  print OUTF "mfba: $MFBA\n";
  print OUTF "rlambda: $RLAMBDA\n";
  print OUTF "alambda: $ALAMBDA\n";
  if ($sieveType==1) {
    print OUTF "q0: $q0\n";
    print OUTF "qintsize: $qrSize\n"; 
    print OUTF "#q1:$q1\n";
  } else {
    print OUTF "a0: $A0\n";
    print OUTF "a1: $A1\n";
    print OUTF "b0: $B0\n";
    print OUTF "b1: $B1\n";
  }
  close(OUTF); 
}

######################################################
sub readParams {
######################################################

  # Read default parameters first. Then we will just override
  # by any user-supplied parameters.
  open(PF, $NAME.".poly");
  my @thisData=<PF>;
  close(PF);
  my @NLINE=grep(/^n:/, @thisData);
  $NLINE[0] =~ /n:\s*(\d+)/;
  $N=$1;

  # Find the polynomial degree.
  my %COEFVALS = ();
  my @COEFLINE=grep(/^c\d+:/, @thisData);
  $D=0;
  while($_ = shift @COEFLINE) {
    # Grab the coefficient index
    # First char of line 'c' followed by a digit string.
    my ($key, $val) = /^(c\d+):\s*(-?\d+)/;
    $key =~ s/c//;
    if ($key > $D) { $D=$key; }
    $COEFVALS{$key} = $val;
  }
  $DEGREE=$D;
#  foreach my $key (reverse sort keys %COEFVALS) {
#    print "-> c$key: $COEFVALS{$key}\n";
#  }

  my @TYPELINE=grep(/type:/, @thisData);
  $TYPE=$TYPELINE[0];
  $TYPE =~ s/.*type: //;
  chomp $TYPE; $TYPE =~ s/\r+$//;
  if ($TYPE =~ /snfs/) {
    # We need the difficulty level of the number, which may be
    # noticably larger than the number of digits.
    my $polyval = new Math::BigInt '0';
    my @MLINE=grep(/^m:/, @thisData);
    $MLINE[0] =~ /m:\s*(\d+)/;
    $M=$1;
    unless(defined $M) {
      @MLINE=grep(/^Y1:/, @thisData);
      $MLINE[0] =~ /Y1:\s*(-?\d+)/;
      my $denom = new Math::BigInt $1;
      @MLINE=grep(/^Y0:/, @thisData);
      $MLINE[0] =~ /Y0:\s*(-?\d+)/;
      my $numer = new Math::BigInt $1;
      $numer = -$numer; 
#      print "-> Common root is $numer / $denom\n";
      my $subtotal = new Math::BigInt;
      for my $i (0..$DEGREE) {
        $subtotal = $COEFVALS{$i} * $numer**$i * $denom**($DEGREE - $i);
        $polyval->badd($subtotal);
      }
    } else {
#      print "-> Common root is $M\n";
      for my $i (reverse 0..$DEGREE) {
        $polyval->bmul($M);
        $polyval->badd($COEFVALS{$i});
      }
    }
#    print "-> Evaluated value is $polyval\n";
    $SNFS_DIFFICULTY = (new Math::BigFloat $polyval)->blog(10);

    printf "-> SNFS_DIFFICULTY is about $SNFS_DIFFICULTY.\n";
    loadDefaultParams($SNFS_DIFFICULTY->bstr(), "snfs");
  } elsif ($TYPE =~ /gnfs/) {
    loadDefaultParams(length($N), $TYPE);
  } else {
    printf "-> Error: poly file should contain one of the following lines:\n";
    printf "-> type: snfs\n";
    printf "-> type: gnfs\n";
    printf "-> Please add the appropriate line and re-run.\n";
    exit;
  }

  # Now look for user-supplied parameters.
  $Q0=0;
  open(PF, $NAME.".poly");
  while (<PF>) {
    chomp;
    s/#.*//; # Remove comments
    s/\s*//g; # Remove whitespace.
    @_ = split /:/;
    my $token=$_[0]; my $val=$_[1];
    if ((length($token)>0) && (length($val)>0)) {
      if ($token eq "n") { $N=$val; }
      elsif ($token eq "m") { $M=$val; }
      elsif ($token eq "rlim") { $RLIM=$val; }
      elsif ($token eq "alim") { $ALIM=$val; }
      elsif ($token eq "lpbr") { $LPBR=$val; }
      elsif ($token eq "lpba") { $LPBA=$val; }
      elsif ($token eq "mfbr") { $MFBR=$val; }
      elsif ($token eq "mfba") { $MFBA=$val; }
      elsif ($token eq "rlambda") { $RLAMBDA=$val; }
      elsif ($token eq "alambda") { $ALAMBDA=$val; }
      elsif ($token eq "knowndiv") { $KNOWNDIV=$val; }
      elsif ($token eq "skew") { $SKEW=$val; }
      elsif ($token eq "q0") { $Q0=$val; }
      elsif ($token eq "qintsize") { $QSTEP=$val; }
      elsif (($token =~ /c./) || ($token =~ /Y./)) {
        push(@COEF, $token.":");
        push(@COEF, $val);
        $COEFHASH{$token} = $val;
      }
    }
  }
  if ($KNOWNDIV) { $KNOWNDIV="-knowndiv ".$KNOWNDIV; }
  if ($Q0==0) {
    $Q0 = $ALIM/2;
  }
  $QSTART=$Q0;
  checkParamFile;
}

######################################################
sub setup {
######################################################

  # Should we resume from an earlier run? #
  $resume=-1;
  if ((-e $RELSBIN.".0")||(-e $JOBNAME)) {
    if ($PROMPTS) {
      print "-> It appears that an earlier attempt was interrupted in progress. Resume? (y/n) ";
      do {
        $_=getc;
        if (($_ eq "Y") || ($_ eq "y")) { $resume=1; }
        elsif (($_ eq "N") || ($_ eq "n")) { $resume=0; }
      } while ($resume < 0);
      printf "\n";
    } else { $resume=1; }
  }

  ############################################
  # Setup the parameters for sieving ranges. #
  ############################################
  if ($resume != 1) {
    if ($CLIENT_ID == 1) {
      # Clean up any junk leftover from an earlier attempt.
      unlink $LOGFILE, <cols*>, <deps*>, 'factor.easy', <lpindex*>;
      unlink <rels*>, 'spairs.out', 'spairs.save.gz', $NAME.'.fb', <*.afb.0>;
      unlink '.params';

      # Was a discriminant divisor supplied?
      if ($DD) {
        open(INFO, $LOGFILE);
        print INFO "$DD\n";
        close(INFO);
      }
      # Create the factor base.
      $cmd="$NICE \"$MAKEFB\" -rl $RLIM -al $ALIM -lpbr $LPBR -lpba $LPBA $LARGEPRIMES -of $NAME.fb -if $NAME.poly";
      print "=>$cmd\n" if($ECHO_CMDLINE);
      $res=system($cmd);
      die "Return value $res. Terminating...\n" if ($res);

    }
    $Q0=$QSTART;
  } else {
    # Get the Q0 value from tmp.job and just restart from there.
    if (-e $JOBNAME) {
      open(INFO, $JOBNAME);
      @lastJobF=<INFO>;
      close(INFO);
      @Q0LINE=grep(/q0:/, @lastJobF);
      $Q0=$Q0LINE[0];
      $Q0 =~ s/.*q0://;
      chomp $Q0;
      $Q0 =~ s/\s//g;
      if (!($Q0)) {
        printf "=> Could not recover q0: field from job file $JOBNAME!\n";
        $Q0=$QSTART;
      }
    } else {
      print "-> File $JOBNAME does not exist. Could not determine a starting q value!\n";
      print "-> Please enter a starting point for the special q: ";
      open(INFO, '-');
      $Q0=<INFO>;
      close(INFO);
    }
  }
}

###############################################################
sub classicalSieve {
# arg0 = A value, to sieve [-A,A]
# arg1 = B0
# arg2 = B1, to sieve for b values in [B0, B1].

  if (!($DOCLASSICAL)) { return };


  if ($#_ < 2) { die("classicalSieve() : Not enough arguments!\n") };
  my $A   = $_[0];
  my $B0  = $_[1];
  my $B1  = $_[2];
  my $maxB = 0;

  # First, scan the log file and find the largest b-value that was sieved.
  open(LOG, $LOGFILE) || return;
  while ($_ = <LOG>) {
    chomp;
    tr/a-z/A-Z/; # Convert to upper case.
    if ( /CLASSICAL/) {
      s/\[.*]\s*CLASSICAL SIEVED//g;  # Strip out everything but the range
      s/\s//g; # Remove whitespace.
      s/\[.*]X//g;  # Strip out the A-range
      s/\[//g; 
      s/]//g; 
      @_ = split /,/;
      if ($_[1] > $maxB) { $maxB = $_[1];}
    }
  }
  $maxB = $maxB + 1;
  if ($maxB < $B0) { $maxB = $B0; }
  if ($maxB >= $B1) { return ;} 
  printf "-> Log file scanned: resuming classical sieve from b=$maxB.\n";
  
  makeJobFile($JOBNAME, 0, 0, 0, 1, -$A, $A, $maxB, $B1);
  $cmd="$NICE \"$CLSIEVE\" -fb $NAME.fb -j $JOBNAME";
  print "=>$cmd\n" if($ECHO_CMDLINE);
  $res=system($cmd);
  unlink $JOBNAME;
  die("Interrupted. Terminating...\n") if ($res);
}



############################################
######## Begin execution here ##############
############################################

$SIG{'INT'}=\&sigDie;
print 
"-> ___________________________________________________________
-> |        This is the factLat.pl script for GGNFS.          |
-> | This program is copyright 2004, Chris Monico, and subject|
-> | to the terms of the GNU General Public License version 2.|
-> |__________________________________________________________|\n";

if (($#ARGV != 0) && ($#ARGV != 2)) {
  print "USAGE: $0 <polynomial file | number file> [ id  num]\n";
  print "  where <polynomial file> is a file with the poly info to use\n";
  print "  or <number file> is a file containing the number to factor.\n";
  print "  Optional: id/num specifies that this is client <id> of <num>\n";
  print "            clients total (clients are numbered 1,2,3,...).\n";
  print " (Multi-client mode is still very experimental - it should only\n";
  print "  be used for testing, and is only intended as a hack for running\n";
  print "  conveniently on a very small number of machines - perhaps 5 - 10).\n";
  exit;
}
$NAME=$ARGV[0];
if ($#ARGV == 2) {
  $CLIENT_ID=$ARGV[1];
  $NUM_CLIENTS=$ARGV[2];
  if (($CLIENT_ID < 1) || ($CLIENT_ID > $NUM_CLIENTS)) {
    print "-> Error: client id should be between 1 and the number of clients ($NUM_CLIENTS)\n";
    exit -1;
  }
} else {
  # Single processor.
  $NUM_CLIENTS=1;
  $CLIENT_ID=1;
}

printf "-> This is client $CLIENT_ID of $NUM_CLIENTS\n";

$NAME =~ s/\.poly//;
$NAME =~ s/\.n//;
print "-> Working with NAME=$NAME...\n";
$JOBNAME=$NAME.".job";
$SIEVER_OUTPUTNAME="spairs.out";
if ($CLIENT_ID > 1) {
  $JOBNAME .= ".$CLIENT_ID";
  $SIEVER_OUTPUTNAME .= "$CLIENT_ID";
}
$psTime=0;

# Is there a poly file already, or do we need to create one?
if (!(-e $NAME.".poly")) {
  printf("-> Error: Polynomial file $NAME.poly does not exist!\n");
  if ($NUM_CLIENTS > 1) {
    printf "-> Script does not support polynomial search across multiple clients!\n";
    exit ;
  }
  unlink $PARAMFILE;
  if (-e $NAME.".n") {
    open(IF, "$NAME".".n");
    my @numberinf=<IF>;
    close(IF);
    my @TMP=grep(/^n:/, @numberinf);
    $TMP[0] =~ /n:\s*(\d+)/;
    $N=$1;
    if (length($N)>0) {
      printf("-> Found n=$N.\n");
      printf("-> Attempting to run polyselect...\n");
      if (length($N)<100) { $USE_KLEINJUNG_FRANKE_PS=0; }
      $psTime=time;
      if ($USE_KLEINJUNG_FRANKE_PS) {
        runPol5;
      } else {
        runPolyselect;
      }
      $psTime=(time-$psTime)/3600;
      die("Polynomial selection failed.\n") unless (-e $NAME.".poly");
    } else {
      printf("-> Could not find a number in the file $NAME.n\n");
      printf("-> Did you forget the 'n:' tag?\n");
      exit;
    }
  }
}

# Read and verify parameters for the factorization.
if (!(-e $NAME.".poly")) {
  die("Cannot find $NAME.poly\n");
}
readParams;
checkParams;
if (!(-e $DEPFILE)) {  setup; }

######################################################################
## Find the minimum number of FF's to get.                          ##
## For this, we must determine the RFB,AFB sizes from the log file. ##
######################################################################
open(INFO, $LOGFILE);
while (<INFO>) {
  @_ = split;
  if ($_[2] eq "RFBsize:") { $rfbSize = $_[3]; }
  if ($_[2] eq "AFBsize:") { $afbSize = $_[3]; }
}
close(INFO);
$totalSize=$rfbSize+$afbSize+64;
$minFF=int($totalSize*(1.0 + $minextraFF));
print "-> minimum number of FF's: $minFF\n";
open(OF, ">>$LOGFILE");
print OF "->               minimum number of FF's: $minFF\n";
close(OF);


# Do some classical sieving, if needed/applicable.
# This is broken - it is still a work in progress!
 classicalSieve($classicalA, 1, $classicalB);
 rename "spairs.out","spairs.add"; 

####################################################
# Finally, sieve until `matbuild' creates a        #
# spmat file, signalling that sieving is done and  #
# we have reached the desired min # of FF's.       #
####################################################
while (!(-e $SPMAT)) {
  printf "-> Q0=$Q0, QSTEP=$QSTEP.\n";
  # Create a job file.
  makeJobFile($JOBNAME, $Q0, $QSTEP, $CLIENT_ID, $NUM_CLIENTS);
  printf("-> Lattice sieving q-values from q=$Q0 to %d.\n", $Q0+$thisQRSize);
  if ($Q0 >= 2**$LPBA) {
    printf "-> $0 : Severe error!\n";
    printf("->     Current special q=$Q0 has exceeded max. large alg. prime = %d !\n", 2**$LPBA);
    printf("-> You can try increasing LPBA, re-launch this script and cross your fingers.\n");
    printf("-> But be aware that if you're seeing this, your factorization is taking\n");
    printf("-> much longer than it would have with better parameters.\n");
    exit;
  }
  $startTime = time;

  # It's very important to call like this, so that if the user CTRL-C's,
  # or otherwise kills the process, we see it and terminate as well.
  $cmd="$NICE \"$LATSIEVER\" -k -o $SIEVER_OUTPUTNAME -v -n$PNUM -a $JOBNAME";
  print "=>$cmd\n" if($ECHO_CMDLINE);
  $res=system($cmd);
  # If we are sieving below the AFB limit, we need to delete the
  # siever's factor base file to make it create a new one. This is
  # a dirty hack - what we should do is modify Franke's code to
  # explicitly allow it.
  if ($sieverAL == ($Q0-1)) {
    unlink <$JOBNAME.afb.*>;
  }
  if ($res) {
    print "-> Return value $res. Updating job file and terminating...\n";
    open(INFO, ".last_spq$PNUM");
    $lastSPQ=<INFO>;
    chomp $lastSPQ;
    close(INFO);
    unlink ".last_spq$PNUM";
    if ($lastSPQ < $Q0) {
      $lastSPQ = $Q0;
    }
    # Move the new relations so they won't be wiped on restart.
    if ($CLIENT_ID == 1) {
      $cmd="\"$CAT\" $SIEVER_OUTPUTNAME >> spairs.add";
      print "=>$cmd\n" if($ECHO_CMDLINE);
      system($cmd);
    } else {
      $cmd="\"$CAT\" $SIEVER_OUTPUTNAME >> spairs.add.$CLIENT_ID";
      print "=>$cmd\n" if($ECHO_CMDLINE);
      system($cmd);
    }
    unlink "$SIEVER_OUTPUTNAME";
    # And update the job file accordingly:
    makeJobFile($JOBNAME, $lastSPQ, $QSTEP, $CLIENT_ID, $NUM_CLIENTS);
    # Record the time to the logfile.
    $stopTime = time;
    $totalTime = $stopTime - $startTime;
    die "Terminating...\n";
  }

  die("Some error ocurred and no relations were found! Examing log file.\n") unless
      -e "$SIEVER_OUTPUTNAME";

  $Q0=$Q1+1;
  if ($CLIENT_ID > 1) {
    $cmd="\"$CAT\" $SIEVER_OUTPUTNAME >> spairs.add.$CLIENT_ID";
    print "=>$cmd\n" if($ECHO_CMDLINE);
    system($cmd);
  } else {
    # Are there relations coming from somewhere else which should be added in?
    if (-e "spairs.add") {
      $cmd="\"$CAT\" spairs.add >> spairs.out";
      print "=>$cmd\n" if($ECHO_CMDLINE);
      system($cmd);
      unlink "spairs.add";
    }
    for ($i=1; $i<=$NUM_CLIENTS; $i++) {
      if (-e "spairs.add.$i") {
        $cmd="\"$CAT\" spairs.add.$i >> spairs.out";
        print "=>$cmd\n" if($ECHO_CMDLINE);
        system($cmd);
        unlink "spairs.add.$i";
      }
    }
    $stopTime = time;
    $totalTime = $stopTime - $startTime;
    $cmd="$NICE \"$PROCRELS\" -fb $NAME.fb -prel $RELSBIN -newrel spairs.out";
    print "=>$cmd\n" if($ECHO_CMDLINE);
    $res=system($cmd);
    die "Return value $res. Terminating...\n" if ($res);
    if ($SAVEPAIRS) {
      $cmd="$NICE \"$GZIP\" -c spairs.out >> spairs.save.gz"; 
      print "=>$cmd\n" if($ECHO_CMDLINE);
      $res=system($cmd);
      die "Return value $res. Terminating...\n" if ($res);
    }
    unlink "spairs.out";
    ##############################################################################
    # CJM, 3/16/05: New code here. This script now decides whether or not
    # there are enough relations to warrant cycle counting and matrix building.
    ##############################################################################
    # Find out how many Relations and total large primes there are.

# this needs to be fixed.
    open(LOG, $LOGFILE) || return;
    while ($_ = <LOG>) {
      chomp;
      tr/a-z/A-Z/; # Convert to upper case.
      if ( /LARGEPRIMES/) {
        s/\[.*]\s*//;  # Strip out the date.
        s/(LARGEPRIMES: |RELATIONS: )//g;  # Strip out the labels.
        s/\s//g; # Remove whitespace.
        @_ = split /,/;
        $tlp=$_[0];
        $trel=$_[1];
      }
    }
    if (($tlp <= 0) || ($trel <= 0)) {
      print "-> Warning: Failed to read total large primes and/or total rels from log!\n";
    } else {
      print "-> Found $tlp total LP vs. $trel relations.\n";
      if ((($tlp - $trel) < 0.8*$trel)||($FORCECC=="on")) {
        $cmd="$NICE \"$MATBUILD\" -fb $NAME.fb -prel $RELSBIN -maxrelsinff $maxRelsInFF -minff $minFF";
        print "=>$cmd\n" if($ECHO_CMDLINE);
        $res=system($cmd);
        die "Return value $res. Terminating...\n" if ($res);
      }
    }

    plotLP;
    # Find out how many FF's there are.
    open(INFO, $LOGFILE);
    while (<INFO>) {
      chomp;
      if (s/finalFF://) {
        $t=$';
      }
    }
    close(INFO);
    if ($t < $minFF) {
      printf "-> Found $t relation-sets versus minFF=$minFF.\n";
      printf "-> More sieving needed.\n";
      # Remove the `spmat' file so we can do some more sieving.
      # Note: this really shouldnt happen anymore, if all the 
      # command-line args are right.
      unlink $SPMAT;
    }
    unlink $JOBNAME;
  }
}

if ($CLIENT_ID > 1) {
  printf "Client $CLIENT_ID terminating...\n";
  exit 0;
}

###############################
# Obviously, the matrix step. #
###############################
if (!(-e $DEPFILE)) {
  print "-> Doing matrix step...\n";
  $cmd="$NICE \"$MATSOLVE\""." -wt ".$matWtFactor; 
  print "=>$cmd\n" if($ECHO_CMDLINE);
  $res=system($cmd);
  die "Return value $res. Terminating...\n" if ($res);
  die("Some error occurred and matsolve did not record dependencies.\n") unless
     (-e $DEPFILE);
} else {
  printf "-> File 'deps' already exists. Proceeding to sqrt step.\n";
}
#############################################
# Do as many square root jobs as needed to  #
# get the final factorization.              #
#############################################
$depnum=0;
my $done = getPrimes;
while (!($done)) {
  $cmd="$NICE \"$SQRT\" $DISC -fb $NAME.fb -deps $DEPFILE -depnum $depnum $KNOWNDIV";
  print "=>$cmd\n" if($ECHO_CMDLINE);
  $res=system($cmd);
  $done = getPrimes;
  if ($depnum++ >= 32) { $done=1; }
} 

if ($CLEANUP) {
  unlink $LOGFILE, <cols*>, <deps*>, 'factor.easy', <lpindex*>;
  unlink <rels*>, 'spairs.out', 'spairs.save.gz', $NAME.'.fb', <*.afb.0>;
  unlink "tmpdata.000";
  unlink $PARAMFILE;
}

# Figure the time scale for this machine.
printf "-> Computing time scale for this machine...\n";
@TMP=`"$PROCRELS" -speedtest`;
@TMP=grep(/timeunit:/,@TMP);
$TIMESCALE=$TMP[0];
$TIMESCALE =~ s/timeunit: //;

# And gather up some stats.
$sieveT=0.0;
$relprocT=0.0;
open(INFO,$LOGFILE);
while (<INFO>) {
  chomp;
  @_ = split;
  if ($_[0] eq "LatSieveTime:") { $sieveT += $_[1]; }
  if ($_[2] eq "RelProcTime:") { $relprocT += $_[3]; }
  if ($_[2] eq "BLanczosTime:") { $matT = $_[3]; }
  if ($_[2] eq "sqrtTime:") { $sqrtT = $_[3]; }
  if (/rels:/) { $rels=$_[2]." ".$_[4]; }
  if (/Initial matrix/) { s/\[.*\] Initial matrix is //; $initmat=$_; }
  if (/Matrix pruned/) { s/\[.*\] Matrix pruned to //; $prunedmat=$_; }
  if (/GGNFS-/) { $version=$_[2]; }
  if (/p: /) {
    s/.*p: //;
    if ((length($_) > 1) && (length($_) < length($N))) {
      push(@DIVISORS, $_);
    }
  }

}
close(INFO);
# Convert times from seconds to hours.
$sieveT /= 3600.0; $relprocT /= 3600.0;
$matT /= 3600.0;   $sqrtT /= 3600.0;
$totalT=$sieveT+$relprocT+$matT+$sqrtT+$psTime;

open($std, ">&STDOUT");         # Save STDOUT handle

if ($TYPE =~ /gnfs/) {
  my $tmpL = length($N);
  my ($dname, $lname) = ($NAME =~ m|(.*/)?(.+)|);
  $sumName="${dname}g${tmpL}-${lname}.txt";
} else {
  my $tmpL=$SNFS_DIFFICULTY->bfloor();
  my ($dname, $lname) = ($NAME =~ m|(.*/)?(.+)|);
  $sumName="${dname}s${tmpL}-${lname}.txt";
}

printf "sumName = $sumName\n";
open(STDOUT, ">$sumName");
print "Number: $NAME\n";
print "N=$N\n";
printf("  ( %d digits)\n", length($N));
if ($TYPE =~ /snfs/) { 
  printf("SNFS difficulty: %d digits.\n", $SNFS_DIFFICULTY);
}
print "Divisors found:\n";

# Sort ascending numerically
@DIVISORS = sort {$a <=> $b} (@DIVISORS);
$r = 1;
while($_ = shift @DIVISORS) {
printf(" r%d=%s\n", $r++, $_);
}

printf("Version: $version\n");
printf("Total time: %1.2f hours.\n", $totalT);
printf("Scaled time: %1.2f units (timescale=%1.3lf).\n", $totalT*$TIMESCALE,$TIMESCALE);
print "Factorization parameters were as follows:\n";
open(PARS, "$NAME.poly");
while (<PARS>) { print $_; }
close(PARS);
print "Factor base limits: $RLIM/$ALIM\n";
print "Large primes per side: $LARGEP\n";
print "Large prime bits: $LPBR/$LPBA\n";
print "Sieved special-q in [$QSTART, $Q0)\n";
print "Relations: $rels\n";
print "Initial matrix: $initmat\n";
print "Pruned matrix : $prunedmat\n";
if ($psTime > 0) {
  printf("Polynomial selection time: %1.2f hours.\n", $psTime);
}
printf("Total sieving time: %1.2f hours.\n", $sieveT);
printf("Total relation processing time: %1.2f hours.\n", $relprocT);
printf("Matrix solve time: %1.2f hours.\n", $matT);
printf("Time per square root: %1.2f hours.\n", $sqrtT);

if ($TYPE =~ /snfs/) { 
  $DIGS = $SNFS_DIFFICULTY->bfloor(); 
  $DEFLINE="$TYPE,$DIGS,$DEGREE,0,0,0,0,0,0,0,0,$RLIM,$ALIM,$LPBR,$LPBA,$MFBR,$MFBA,$RLAMBDA,$ALAMBDA,$QINTSIZE\n";
}
else {
  $DIGS=int(log($N)/log(10.0));
  $DEFLINE="$TYPE,$DIGS,$DEGREE,maxs1,maxskew,goodScore,efrac,j0,j1,eStepSize,maxTime,$RLIM,$ALIM,$LPBR,$LPBA,$MFBR,$MFBA,$RLAMBDA,$ALAMBDA,$QINTSIZE\n";
}
printf("Prototype def-par.txt line would be:\n$DEFLINE");
printf("total time: %1.2f hours.\n", $totalT);
print " --------- CPU info (if available) ----------\n";
close(STDOUT);


if (-x '/bin/dmesg') {
  $cmd="/bin/dmesg | grep CPU | grep stepping >> $sumName";
  print "=>$cmd\n" if($ECHO_CMDLINE);
  system($cmd);
  $cmd="/bin/dmesg | grep Memory >> $sumName";
  print "=>$cmd\n" if($ECHO_CMDLINE);
  system($cmd);
  $cmd="/bin/dmesg | grep -i bogomips >> $sumName";
  print "=>$cmd\n" if($ECHO_CMDLINE);
  system($cmd);
}
if (-x '/usr/sbin/x86info') {
  $cmd="/usr/sbin/x86info -mhz >> $sumName";
  print "=>$cmd\n" if($ECHO_CMDLINE);
  system($cmd);
}

open(STDOUT, ">&", $std); # Restore orig STDOUT
print "-> Factorization summary written to $sumName.\n";

# Send 5 BELLS and flush output so \n not necessary (in Cygwin) after each BELL
use IO::Handle;
STDOUT->autoflush(1);

for($i = 0; $i < 5; $i++){sleep 1;print "\a";}