? count=0;for(W=1,5555130,for(Z=W,W, par=W;new=3; for(X=1,par,new=4*new+1);if( ispseudoprime(new,1) ,count=count+1; print(W,” “,count,” “,round(log(new)/log(10)) ) )) )
1 1 1
2 2 2
4 3 3
5 4 4
8 5 5
10 6 7
11 7 7
20 8 13
32 9 20
40 10 25
131 11 79
257 12 155
263 13 159
350 14 211
448 15 270
634 16 382
725 17 437
803 18 484
832 19 501
1769 20 1066
2363 21 1423
2548 22 1535
3544 23 2134
4508 24 2715
4904 25 2953
5494 26 3308
N.B.: strong Rabin-Miller test to one base only, for more speed.
