- --- In primenumbers@yahoogroups.com,

"mikeoakes2" <mikeoakes2@...> wrote:

> Any pari code I come up with is horribly sluggish

This was the puzzle:

> for n more than about 10^6.

> I await your 4 lines with interest!

> Puzzle: Find a composite positive odd integer n that passes

n = 7056721 is the unique solution with n < 10^10.

> the Lucas test V(a,1,n) = a mod n, for every integer a.

However, I have more than 144000 larger solutions on file.

Using Pari-GP, I found a candidate solution as follows:

tn(n,a)=2*polcoeff(lift(Mod((a+x)/Mod(2,n),x^2+4-a^2)^n),0)==Mod(a,n);

ts(n)=tn(n,3)&&tn(n,4)&&tn(n,5)&&tn(n,6)&&tn(n,7)&&tn(n,8);

{forstep(n=3,10^8,2,if(n%3&&!isprime(n)&&ts(n),

print(n);print(round(gettime/10^3)" seconds");break()));}

7056721

84 seconds

Then by googling> pseudoprime 7056721

I arrived at

http://adm.lnpu.edu.ua/downloads/issues/2008/N2/adm-n2-4.pdf

which proves that n is a solution if and only if it is an odd

square-free composite integer such that for each prime p|n

n = +/- 1 mod p-1 ... [1]

n = +/- 1 mod p+1 ... [2]

This paper claims that 7056721 is the unique solution with n < 10^10.

It is easy to construct larger solutions. For example

http://tech.groups.yahoo.com/group/primenumbers/message/19971?var=0

provides a link to 246 solutions in

http://physics.open.ac.uk/~dbroadhu/cert/erdos2.out

with n = 1 mod p^2-1, for each prime p|n.

For 144153 such solutions, see the 12 MB file

http://physics.open.ac.uk/~dbroadhu/cert/carmrob2.txt

obtained by Robert Gerbicz on 1 April 2009 and verified here:

read("carmrob2.txt");c=0;

{for(k=1,#ls,n=ls[k];ff=factor(n)[,1];

for(j=1,#ff,p=ff[j];if(n%(p^2-1)!=1,print([n,p]);c++)));}

print(c"/"#ls" failures in "round(gettime/10^3)" seconds");

0/144153 failures in 27 seconds

Puzzle: Find the smallest composite integer n > 7056721

such that V(a,1,n) = a mod n, for every integer a.

Comment: The solution is not known to me.

David Broadhurst, 24 May 2010 - --- In primenumbers@yahoogroups.com, "djbroadhurst" <d.broadhurst@...> wrote:
>

You did very much what I did.

> I tried 1/n^c:

>

> v=[1237.1, 328.7, 105.4, 28.01, 6.22 , 1.510, 0.439, 0.0939];

> print(vector(8,k,-log(v[k]/10^6)/(k+4)/log(10)));

>

> [0.5815, 0.5805, 0.5682, 0.5691, 0.5785, 0.5821, 0.5780, 0.5856]

>

> and then A/n^c, using the first datum to remove A:

>

> print(vector(7,k,-log(v[k+1]/v[1])/k/log(10)));

>

> [0.5756, 0.5348, 0.5484, 0.5747, 0.5827, 0.5750, 0.5885]

>

> In both cases c =~ Euler looks rather convincing,

> given the statistics. Well spotted, Sir!

I nearly fell off the chair when I averaged everything out and saw 0.57... :-)

> How strongly are you committed to A = 1, for the average,

Not very.

> given the variability of the overall factor with a?

Would you buy an appeal to Occam's razor, mon vieux?

Mike