Loading ...
Sorry, an error occurred while loading the content.
 

Re: Lucas super-pseudoprime puzzle

Expand Messages
  • mikeoakes2
    ... The best pari-GP script I can come up with is /quadratic/ in n_max, so will take about 23360000 secs or 9 months at 3.6 GHz to reach n_max=10^7 and
    Message 1 of 46 , Nov 2, 2010
      --- In primenumbers@yahoogroups.com, "djbroadhurst" <d.broadhurst@...> wrote:
      >
      > Definition: The integer pair (q,n) is a "Lucas super-pseudoprime"
      > (LSPS) pair if and only if n is composite, n > q > 0, and
      > V(p,q,n) = p mod n, for every integer p.
      >
      > With n = 11*29*37*43 = 507529, there are precisely 37 LSPS pairs
      > (q,n), namely those with q = (11*k + 3)*29*43, for k = 0 to 36.
      > Mike Oakes remarked that 507529 is not a Carmichael number.
      >
      > Puzzle: Find an odd square-free non-Carmichael number
      > such that there exist more than 70,000 LSPS pairs (q,n).
      >
      > Comment: There is a solution with n < 10^7.

      The best pari-GP script I can come up with is /quadratic/ in n_max, so will take about 23360000 secs or 9 months at 3.6 GHz to reach n_max=10^7 and (presumably) solve your puzzle.
      Here it is:-

      \\For odd non-prime n, finds first q s.t. lucasV(p,q,n)=p mod n for all p
      \\(Using CRT, q can be reconstructed from its remainders)
      v(p,q,n,m)=2*polcoeff(lift(Mod((p+x)/Mod(2,m),x^2+4*q-p^2)^n),0);
      n_min=3; n_max=10^5;
      default(primelimit,n_max);
      gettime;
      {
      forstep(n=n_min,n_max,2,
      if((n%10^5)==1,print("n="n));
      f=factor(n);
      facs=#f[,1];
      n_ok=1;
      if(facs>1,
      for(ind=1,facs,
      if(f[ind,2]>1,n_ok=0;break(1)); \\not squarefree
      m=f[ind,1];
      count=0;
      for(q=1,m,
      ok=1;
      for(p=1,m,
      if(lift(v(p,q,n,m)-p)<>0,ok=0;break(1));
      ); \\end for p
      if(ok,count++;break(1)); \\add code here to remember q, for CRT
      ); \\end for q
      if(count==0,n_ok=0;break(1));
      ); \\end for ind
      if(n_ok,print("n="n" "factor(n)" is ok"));
      ); \\end if facs>1
      ); \\end forstep n
      tim=gettime;
      print(" took "tim" ms");
      }
      /*
      n_max time
      10^4 24 secs
      10^5 2336 secs
      */

      What am I missing??

      Mike
    • djbroadhurst
      ... http://physics.open.ac.uk/~dbroadhu/cert/dbmo116.out gives my 116, in the format [n, factors, number of solutions] With n
      Message 46 of 46 , Nov 9, 2010
        --- In primenumbers@yahoogroups.com,
        "mikeoakes2" <mikeoakes2@...> wrote:

        > > My revised count up to 2*10^10 is 116.
        > My (original) count up to 2*10^10 was 105.
        > So it must have missed 11, i.e. a bigger proportion.

        http://physics.open.ac.uk/~dbroadhu/cert/dbmo116.out
        gives my 116, in the format [n, factors, number of solutions]

        With n < 2*10^10, the record-holder for the number of solutions is
        [2214495361, [13, 17, 23, 29, 83, 181], 147407]
        which googles quite nicely, linking to
        http://www.cs.rit.edu/usr/local/pub/pga/fibonacci_pp

        David
      Your message has been successfully submitted and would be delivered to recipients shortly.