|
PROGRAM
|
(PARI) \\ uses Pari functions ispseudoprime, denominator, and list \\ tp2cnt =count of adjacent, unique twin primes (lower) \\ p2cnt = count of unique adjacent primes \\ N = Bernoulli # maxx= max N idx*=list pointers \\ u=array of adj. unique TP u2=array of adj. unique primes cnt=0; pcnt=0; n=1; opcnt=0; plcnt=0; phcnt=0; tpcnt=0; otpcnt=0; tp2cnt=0; p2cnt=0; limmit=429000001; \\4E8 idx=1; \\init idx2=1;
\\init if(maxx<2, maxx=2); if(maxx>200000, maxx=200000); print ("enforced n max = ", maxx); u=listcreate(500); \\unique entries < 10% TP u2=listcreate(5000); \\unique entries < 50% while(n<maxx+1, if(bernfrac(n)==0, n++) ; \\guard against 0 division m=abs(denominator(bernfrac(n)) ); p=m-1; q=m+1;
\\ we check to either side cnt++; if( ispseudoprime(p)==1, plcnt++ ); if( ispseudoprime(q)==1, phcnt++ ); if( ispseudoprime(q)==1 || ispseudoprime(p)==1, pcnt++ ); if( ispseudoprime(q)==1 && ispseudoprime(p)==1, tpcnt++ ); if( tpcnt>otpcnt, print("TP! @ ", n, " ", q, " "
\\IF.FOR, IF if(p>0, listput(u, p, idx); tp2cnt++ ); if(p>0, idx++) ); \\ end IF if( pcnt>opcnt, if(p>limmit, p=0); \\lower if(ispseudoprime(p)==0, p=0); if(p>0, for(z=1, idx2-1, if(p==u2[z], p=0); ); ); if(p>0, listput(u2, p, idx2); p2cnt++); if(p>0, idx2++); if(q>limmit, q=0);
\\upper if(ispseudoprime(q)==0, q=0); if(q>0, for(z=1, idx2-1, if(q==u2[z], q=0); ); ); if(q>0, listput(u2, q, idx2); p2cnt++); if(q>0, idx2++) ); \\ end IF opcnt=pcnt; otpcnt=tpcnt; n++ ); \\end WHILE if(cnt>0, print(2.*p2cnt/cnt) ); \\reflect every other (n.z.) print ("TP count = ", tp2cnt ); print ("enforced n max = ", maxx);
|
