Possible bug in __gmpn_tdiv_qr

[Michael Reilly]

-----BEGIN PGP SIGNED MESSAGE-----
Hash: SHA1

Hi there, I reached a backtrace while running a program, and I was
wondering if I had uncovered a bug.

There are really three things this could be:
1) User error
2) Computer error (RAM error maybe)
3) Library error

So I figured I'd do a dump here and let you know what I came across,
and let you decide if it's my or my computer's fault.

First the bt:

*** glibc detected *** ./arithprogalgorithmic: double free or
corruption (!prev): 0x0000000001ece700
***                                                                      
 
 
======= Backtrace:
=========                                                                                                                                                
 
 
 
/lib/libc.so.6[0x7fcc601d8cb8]
/lib/libc.so.6(cfree+0x76)[0x7fcc601db276]
/usr/lib/libgmp.so.3(__gmp_tmp_reentrant_free+0x22)[0x7fcc5ff27432]  
/usr/lib/libgmp.so.3(__gmpn_tdiv_qr+0x407)[0x7fcc5ff51ec7]
/usr/lib/libgmp.so.3(__gmpq_get_d+0x12c)[0x7fcc5ff41c5c]
./arithprogalgorithmic[0x40232b]  
/lib/libc.so.6(__libc_start_main+0xe6)[0x7fcc6017f5a6]
./arithprogalgorithmic[0x401bc9]

Next the code (compiled with -O3):

#include <iostream>
#include <iomanip>
#include <fstream>
#include <cmath>
#include <gmpxx.h>
#include <string>
#include <sstream>
#define MAXITER 10000000
#define MAXPRIME 1000  

int main() {
  long int primes[MAXPRIME];
  primes[0] = 2;           
  long int numfound = 1;   
  long int nextprime = 3;  
  for (int i=1; i<MAXPRIME; i++) {
    for (int j=0; j<numfound; j++) {
      if (nextprime % primes[j] == 0) {
        j=0;                          
        nextprime +=2;                
      }                               
    }                                 
    primes[numfound]=nextprime;       
    numfound++;                       
  }                                   

  for (long int p=1; p<MAXPRIME; p++) {
    mpq_class sum;                    
    sum = 1;                          
    long int series_length = 1;       
    long int length_arith_prog = primes[p]-1; //only valid for numbers
one less than an odd prime
    long int n = 1;
    std::cout << std::setprecision(20);
    std::ofstream results;
    std::ostringstream ostr;
    std::string filename;
    ostr << "results-k" << length_arith_prog << "-algo.txt";
    results.open(ostr.str().c_str());
    results << std::setprecision(20);

    for (;series_length<MAXITER;) {
      long int powerofn = 0;
      long int temp = series_length;
      while (temp%length_arith_prog == 0) {
        temp /= length_arith_prog;
        powerofn++;
      }
      long int diff = 1;
      for (long int i=0; i<powerofn; i++) {
        diff += std::pow(2*length_arith_prog-1,i);
      }
      n += diff;
      series_length++;
      mpq_class tempsum(1,n);
      sum += tempsum;
      results << n << " " << series_length << " " << std::setw(5) <<
diff << " ";
      results << sum.get_d() << std::endl;
    }
    results.close();
    std::cout << "Finished k=" << length_arith_prog << std::endl;
  }
}

Next the output (program died about 6 days into calculation on my
AMD64 X2) the full files are available if needed:
% tail results-k4-algo.txt
301152635 3536109     2 5.6628709140244914977
301152636 3536110     1 5.662870917345067312
301152637 3536111     1 5.6628709206656422381
301152638 3536112     1 5.6628709239862171643
301152647 3536113     9 5.6628709273067929786
301152648 3536114     1 5.6628709306273679047
301152649 3536115     1 5.6628709339479428309
301152650 3536116     1 5.662870937268517757
301152652 3536117     2 5.6628709405890926831
301152653 3536118     1 5.6628709439096684974

Finally some information about my system and the version of libgmp used.
gcc version 4.3.3-5ubuntu4
libgmp 4.3.1 (compiled with make check, self reported working.)
Linux oberus 2.6.28-14-generic #46-Ubuntu SMP Wed Jul 8 07:41:18 UTC
2009 x86_64 GNU/Linux
./config.guess = athlon64-unknown-linux-gnu
./configfsf.guess = x86_64-unknown-linux-gnu
-----BEGIN PGP SIGNATURE-----
Version: GnuPG v1.4.9 (GNU/Linux)
Comment: Using GnuPG with Mozilla - http://enigmail.mozdev.org

iQEcBAEBAgAGBQJKdnujAAoJEG2iyL6HyDY6kYAH/ifG2i07Q0MFfVcrdn/exuPJ
auM23hQp06kuTdclT+5+mg5W0NzHTJ7MWQhTWdErZ/viEN5wM+P6KzQBfY/Gb9uI
ORwuDQb9xUVbiE8dq1FEsnEmVLomxMnO5b75KmzsDcGkl2syVJOjxSMBrZ1W8RKt
Abth6ps6QMHFfB9+rd5IQFu5wlWLt2DIpq075UsKVGfxIqZrvtTwmYaxhmN6jxW5
JQoA2pbMn2oP2Zbs0psZbZ2COpn+4OO25Kb3BUNk3dIGVO9Yh2holw3besq/xW4d
DcaEq/y2FMGSZhfgqRrRoOWi3dXtLmiC0079dP07yOesP9AMosHUGWGCGf/IA9E=
=/Gk8
-----END PGP SIGNATURE-----