Msieve with GNFS support

[frmky]

The version I use on small projects is basically the version in the GGNFS Yahoo group. (OK, not really. But I've just changed it to wait for more relations before the first msieve run and to support my MPI'ed sievers so I can run it on the cluster.) But I don't stop and restart it. When I have stopped it in the past, it seemed to work OK on restart, but that's not something I've spent time thinking about.

Greg

[FactorEyes]

Originally Posted by frmky

Sat Sep 20 03:10:46 2008  commencing square root phase
Sat Sep 20 03:10:46 2008  reading relations for dependency 1
Sat Sep 20 03:11:00 2008  read 8759041 cycles
Sat Sep 20 03:16:10 2008  cycles contain 28802530 unique relations
Sat Sep 20 03:25:11 2008  read 28802530 relations
Sat Sep 20 03:54:55 2008  multiplying 41018894 relations
Sat Sep 20 09:49:26 2008  multiply complete, coefficients have about 1183.37 million bits
Sat Sep 20 09:50:10 2008  initial square root is modulo 203383
Sun Sep 21 00:46:01 2008  prp76 factor: 2975299737951435964540299996898411217993568426081224805202105572310472987267
Sun Sep 21 00:46:01 2008  prp189 factor: 137741818854348583884185387895829731532755157992504981233434514906522897328661266349593048714205323208173427503632881246417304040292597529414966939921057660532055795650129085411468978166003
Sun Sep 21 00:46:01 2008  elapsed time 656:26:13

/unhelpful/
Hmmm. If the splitting dependency were 4 or 5, this would be more typical of my runs. At 25 hours per square root, that's processor-tastic! I think I need the secret of making that first dependency count.
/\unhelpful/

[miklin]

Originally Posted by Batalov

One thing I noticed (for the iterative reruns of the matrix stage) is that the *.fb file gets slightly bigger after an attempted matrix build and later cannot be used (or the *.dat file will be emptied by the binary -- something somewhere doesn't match and that signals msieve that this is a new project, so off everything goes). Free relations? But these are probably added to *.dat file, not to *.fb. This may be caused by the factMsieve.pl script and not the binary. Is this a known effect?

--Serge

If that is necessary on what for the reasons to stop the program, stop. Then msieve.dat it is necessary to cut out from a file the first line with number (N ).
Further it is necessary to rename the received file msieve.dat--> spairs.add.
That's all.
If after restart at you msieve.dat will change the size there where at you the file spairs.out has started to be created.
There you add spairs.add. Also do not endure further a script will make all.

I so do. English at me bad excuse.

Sergey Miklin

#include <stdio.h>

//------------------------------------------------------------------------------------------------------
void remoteN(char* filein, char* fileout, char* filen)
{
FILE* fin;                                                                                                
fin=fopen(filein, "rt");
if(fin==NULL) {fprintf(stdout, "can not open %s\n", filein); return;}

FILE* fout;
fout=fopen(fileout, "wt");
if(fout==NULL) {fprintf(stdout, "can not open %s\n", fout); fclose(fin); return;}

FILE* fn;
fn=fopen(filen, "wt");
if(fn==NULL) {fprintf(stdout, "can not open %s\n", filen); fclose(fin); fclose(fout); return;}

// get input file len

fseek(fin, 0, SEEK_END);
size_t FinLen = ftell(fin);
/*
char* WorkDim;
if(FinLen < (256*1024*1024)) WorkDim = new char[FinLen];
else                         WorkDim = new char[256*1024*1024];

if(WorkDim==NULL)
  {
  fprintf(stdout, "can not open memory");
  fclose(fin);
  fclose(fout);
  fclose(fn);
  return;
  }
*/


fseek(fin, 0, SEEK_SET);
char* WorkDim = new char[4096];
while(fgets(WorkDim, 4096, fin)!=NULL)
  {
  if((WorkDim[0] == 'N')||(WorkDim[0]=='n')) 
    {
    if(fputs(WorkDim, fn)==EOF)
      {
      fprintf(stdout, "can not write file %s", filen);
      delete WorkDim;
      fclose(fin);
      fclose(fout);
      fclose(fn);
      return;
      }
    } 
  else if(fputs(WorkDim, fout)==EOF)
    {
    fprintf(stdout, "can not write file %s", fileout);
    delete WorkDim;
    fclose(fin);
    fclose(fout);
    fclose(fn);
    return;
    }
  }

delete WorkDim;
fclose(fin);
fclose(fout);
fclose(fn);
return;
}
//------------------------------------------------------------------------------------------------------
void addN(char* filein, char* fileout, char* filen)
{
FILE* fin;                                                                                                
fin=fopen(filein, "rt");
if(fin==NULL) {fprintf(stdout, "can not open %s\n", filein); return;}

FILE* fout;
fout=fopen(fileout, "wt");
if(fout==NULL) {fprintf(stdout, "can not open %s\n", fout); fclose(fin); return;}

FILE* fn;
fn=fopen(filen, "rt");
if(fn==NULL) {fprintf(stdout, "can not open %s\n", filen); fclose(fin); fclose(fout); return;}

fseek(fin, 0, SEEK_SET);
char* WorkDim = new char[4096];


// Add N from N file
if(fgets(WorkDim, 4096, fn)==NULL)
  {
  fprintf(stdout, "can not read file %s", filen);
  delete WorkDim;
  fclose(fin);
  fclose(fout);
  fclose(fn);
  return;
  }
if(fputs(WorkDim, fout)==EOF)
  {
  fprintf(stdout, "can not write file %s", fileout);
  delete WorkDim;
  fclose(fin);
  fclose(fout);
  fclose(fn);
  return;
  }

while(fgets(WorkDim, 4096, fin)!=NULL)
  {
  if(fputs(WorkDim, fout)==EOF)
    {
    fprintf(stdout, "can not write file %s", fileout);
    delete WorkDim;
    fclose(fin);
    fclose(fout);
    fclose(fn);
    return;
    }
  }
delete WorkDim;
fclose(fin);
fclose(fout);
fclose(fn);
return;
}
//------------------------------------------------------------------------------------------------------
void Massege(void)
{
fprintf(stdout, "remote N from file\n"); 
fprintf(stdout, "format r Nwork <file.in> <file.out> <n.file>\n");
fprintf(stdout, "add N to file\n"); 
fprintf(stdout, "format a Nwork <file.in> <file.out> <n.file>\n");
}
//------------------------------------------------------------------------------------------------------
int main(int argc, char* argv[])
{
if(argc<5)
  {
  Massege();
  return 1;
  }
fprintf(stdout, "\n"); 
 
switch(argv[1][0])
  {
  case 'r':
  case 'R': remoteN(argv[2], argv[3], argv[4]); break;
  case 'a':
  case 'A': addN(argv[2], argv[3], argv[4]); break;
  default: Massege(); return 1;
  }
return 0;
}

[Batalov]

Cracking sqrt, Jason! Kudos!

I've just repeated a small (GNFS-144) sqrt; indeed, almost half the memory (it would have been half, I think, if the odd multiplicity rels were half of all, but there's ~60% of them) and half the time (54min vs. 111minutes).

Just curious: you've now removed the sanity check "square roots are identical". Does it mean that it will not happen anymore? I did see such error once, actually.

Serge

[jasonp]

Originally Posted by Batalov

Cracking sqrt, Jason! Kudos!

[...]

Just curious: you've now removed the sanity check "square roots are identical". Does it mean that it will not happen anymore? I did see such error once, actually.

Thanks for the test; if you want fast, the CADO square root is something like 5x faster than the msieve version, and probably still accelerating.

The 'square roots are identical' error is not an error at all; if you have modular square roots X and Y, and gcd(X+Y, N) is 1 or N, that means X = +-Y mod N. So if X == Y the dependency is just unlucky, nothing is broken.

[R.D. Silverman]

Originally Posted by jasonp

Thanks for the test; if you want fast, the CADO square root is something like 5x faster than the msieve version, and probably still accelerating.

The 'square roots are identical' error is not an error at all; if you have modular square roots X and Y, and gcd(X+Y, N) is 1 or N, that means X = +-Y mod N. So if the dependency doesn't work it's because X and Y are unlucky, not because anything is broken.

And the CWI square root is faster still; I have never seen it take more than
1/2 hour on a dependency.

[xilman]

Originally Posted by R.D. Silverman

And the CWI square root is faster still; I have never seen it take more than 1/2 hour on a dependency.

I've seen it take well over two hours.

Of course, it all depends on, at least, the speed of the machine and the number of relations in the dependency.


Paul

[jasonp]

Originally Posted by xilman

Of course, it all depends on, at least, the speed of the machine and the number of relations in the dependency.

The brute force square root also speeds up when the leading algebraic polynomial coefficient is small, since that reduces convolution sizes.

As long as brute force is in the same ballpark I'm happy. Actually the current version is something of a relief, because it doubles the size of coefficients needed before risking roundoff error problems due to use of floating point FFTs, meaning the code can now handle square root problems twice as large.

[frmky]

Version 1.38 trimmed 11 hours off of the square root, from 21.5 hours to 10.5 hours.

Wed Sep 24 11:34:53 2008  commencing square root phase
Wed Sep 24 11:34:53 2008  reading relations for dependency 1
Wed Sep 24 11:35:12 2008  read 8759041 cycles
Wed Sep 24 11:40:27 2008  cycles contain 28802530 unique relations
Wed Sep 24 11:47:35 2008  read 28802530 relations
Wed Sep 24 12:12:33 2008  multiplying 23401912 relations
Wed Sep 24 15:44:34 2008  multiply complete, coefficients have about 676.91 million bits
Wed Sep 24 15:44:58 2008  initial square root is modulo 1182901
Wed Sep 24 21:58:53 2008  prp76 factor: 2975299737951435964540299996898411217993568426081224805202105572310472987267
Wed Sep 24 21:58:53 2008  prp189 factor: 137741818854348583884185387895829731532755157992504981233434514906522897328661266349593048714205323208173427503632881246417304040292597529414966939921057660532055795650129085411468978166003
Wed Sep 24 21:58:53 2008  elapsed time 10:24:04

[frmky]

Originally Posted by jasonp

It won't cause any problems, but I expect the result to be just about identical except the matrix should end up slightly smaller. The filtering changes in v1.38 will only affect datasets with a large amount of excess relations.

You're exactly right. With the same set of relations, 1.38 produced a matrix that is slightly smaller and lighter.

1.38 produced

Fri Sep 26 04:55:55 2008  matrix is 17160155 x 17160700 (5135.0 MB) with weight 1579620033 (92.05/col)
Fri Sep 26 04:55:55 2008  sparse part has weight 1157330594 (67.44/col)
Fri Sep 26 05:36:36 2008  filtering completed in 3 passes
Fri Sep 26 05:36:52 2008  matrix is 17124041 x 17124241 (5130.7 MB) with weight 1578041823 (92.15/col)
Fri Sep 26 05:36:52 2008  sparse part has weight 1156617685 (67.54/col)
Fri Sep 26 05:43:10 2008  read 17124241 cycles
Fri Sep 26 05:43:56 2008  matrix is 17124041 x 17124241 (5130.7 MB) with weight 1578041823 (92.15/col)
Fri Sep 26 05:43:56 2008  sparse part has weight 1156617685 (67.54/col)
Fri Sep 26 05:43:57 2008  saving the first 48 matrix rows for later
Fri Sep 26 05:44:27 2008  matrix is 17123993 x 17124241 (4915.2 MB) with weight 1209033267 (70.60/col)
Fri Sep 26 05:44:28 2008  sparse part has weight 1117239156 (65.24/col)

while 1.36 produced

Mon Jun 30 00:16:00 2008  matrix is 17551760 x 17552222 (5276.9 MB) with weight 1622323775 (92.43/col)
Mon Jun 30 00:16:00 2008  sparse part has weight 1190223841 (67.81/col)
Mon Jun 30 00:55:27 2008  filtering completed in 3 passes
Mon Jun 30 00:55:39 2008  matrix is 17516334 x 17516534 (5272.8 MB) with weight 1620810005 (92.53/col)
Mon Jun 30 00:55:39 2008  sparse part has weight 1189559600 (67.91/col)
Mon Jun 30 01:01:03 2008  read 17516534 cycles
Mon Jun 30 01:01:33 2008  matrix is 17516334 x 17516534 (5272.8 MB) with weight 1620810005 (92.53/col)
Mon Jun 30 01:01:33 2008  sparse part has weight 1189559600 (67.91/col)
Mon Jun 30 01:01:34 2008  saving the first 48 matrix rows for later
Mon Jun 30 01:02:00 2008  matrix is 17516286 x 17516534 (5050.2 MB) with weight 1242199788 (70.92/col)
Mon Jun 30 01:02:01 2008  sparse part has weight 1148718869 (65.58/col)

Greg

[jasonp]

Cool. Your test case saw a bigger improvement than I expected, probably because your input number generated many more free relations than usual.