PARI's commands

[science_man_88]

Originally Posted by 3.14159

@Sm88: Did you define what the Mersennes are? (Or at least Mersenne_Primes(n)?)

the code for it will figure it out so the best I can do right now(off the top of my head) is:

a=0;for(p=1,x,if(isprime(p) && isprime(2^p-1),a=a+1;print(p);if(a==n,break();))

but I haven't figured out how to do my method suggested earlier involving indexes of A002450 or A121290.

[3.14159]

Originally Posted by CRGreathouse

The naive method would test all ~4e97 bases. Each base would require the calculation of b^(n-1) mod n. These take about a millisecond each, for a total of ~4e94 seconds or ~1e85 lifetimes. Factoring n, if it wasn't written as above, could be done in a few hours. Big win for 'my' method.

Okay. But it would be rendered useless for any number larger than about 90-105 or so digits. I was right about it being unnecessary for larger integers.

Also: @Sm88: Isn't it supposed to be

a=0;forprime(p=1,x,if(isprime(p) && isprime(2^p-1),a=a+1;print(p);if(a==n,break();))

?

(My emphasis in red.)

[science_man_88]

Originally Posted by 3.14159

Also: @Sm88: Isn't it supposed to be

a=0;forprime(p=1,x,if(isprime(p) && isprime(2^p-1),a=a+1;print(p);if(a==n,break();))

?

only if you know before hand x is under primelimit if it isn't forprime may not work as it uses primes[i] last I checked. so if you don't already have all the primes needed calculated it must calculate them on the fly and the isprime(p) would need to be eliminated if you use forprime as it's not needed until after exhausting primelimit.

[CRGreathouse]

Originally Posted by 3.14159

Okay. But it would be rendered useless for any number larger than about 90-105 or so digits. I was right about it being unnecessary for larger integers.

When I said larger integers, I meant integers larger than 2^32, since anything larger will cause problems at word precision for 32-bit computers. If you're using a 64-bit Linux system (the Windows version is 32-bit even on 64-bit computers) then the problem crops up above 2^64 or so. Since those ranges are quite factorable, it *is* relevant.

But even if you were concerned only about numbers below 2^32 ~ 4e9, I still wouldn't recommend using a complicated inelegant method.

[CRGreathouse]

Originally Posted by science_man_88

only if you know before hand x is under primelimit if it isn't forprime may not work as it uses primes[i] last I checked. so if you don't already have all the primes needed calculated it must calculate them on the fly and the isprime(p) would need to be eliminated if you use forprime as it's not needed until after exhausting primelimit.

It uses the list of primes, yes. It doesn't call prime(i) which is actually pretty inefficient 'under the hood' as I recall.

And if you did use forprime, you wouldn't need to have isprime(p), just isprime(2^p-1).

[science_man_88]

Originally Posted by CRGreathouse

It uses the list of primes, yes. It doesn't call prime(i) which is actually pretty inefficient 'under the hood' as I recall.

And if you did use forprime, you wouldn't need to have isprime(p), just isprime(2^p-1).

yes and if x> primelimit it would hand out :not enough primes precomputed or what ever it says. so if x>primelimit the isprime(p) is needed.

this is something the code creating script would have to figure out to give the most efficient codes.

[axn]

An interesting aside:

nu=valuation(f[1]-1,2);
for(i=2,#f,nu=min(nu,valuation(f[i]-1,2))

can be replaced by

nu=valuation(f-vectorv(#f,X,1),2)

vectorv, because f is a column vector (and because PARI doesn't support vector/scalar addition).

[science_man_88]

(12:03) gp > nu=valuation(f-vectorv(#f,X,1),2)
  *** _-_: forbidden addition t_POL + t_COL.

doesn't look like it.

and your first one came back syntax errors then I fixed it and :

(13:10) gp > for(i=2,#f,nu=min(nu,valuation(f[i]-1,2)))
  *** for: _[_]: not a vector.

sorry forgot one thing.

(13:12) gp > nu=valuation(f[1]-1,2);for(i=2,#f,nu=min(nu,valuation(f[i]-1,2)))
  ***   _[_]: not a vector.

[science_man_88]

so

1)read(post)
2)dissect the post into parts we can search addhelp() with for a function that it will solve it
3) post it on the forum ?

[3.14159]

Excellent.

A small inconvenience:

Here is a code snippet that allows me to search for 2nd kind Cunningham semiprimes (Of the form 2p²-p that are pseudoprime to base 2):

modbmpsp(x,n,m)=if(isPRP(modbm(x,n,m),b=2),print(modbm(x,n,m))

After it runs partially: Pari says:

***Mod: division by zero.

Here is my PRP code:

isPRP(n,b=2)={
	my(s=valuation(n-1,2),d=n>>s);
	n=Mod(b,n)^d;
	if(n==1,return(1));
	while(s--,
		if(n==-1,return(1));
		n=n^2
	);
	n==-1
};

Sorry for the question, but: Where is the code asking for a division by zero?

[science_man_88]

Originally Posted by 3.14159

Excellent.

A small inconvenience:

Here is a code snippet that allows me to search for 2nd kind Cunningham semiprimes (Of the form 2p²-p that are pseudoprime to base 2):

modbmpsp(x,n,m)=if(isPRP(modbm(x,n,m),b=2),print(modbm(x,n,m))

After it runs partially: Pari says:

***Mod: division by zero.

Here is my PRP code:

isPRP(n,b=2)={
	my(s=valuation(n-1,2),d=n>>s);
	n=Mod(b,n)^d;
	if(n==1,return(1));
	while(s--,
		if(n==-1,return(1));
		n=n^2
	);
	n==-1
};

Sorry for the question, but: Where is the code asking for a division by zero?

my guess is since n is not initiated, it is taken as 0 in which case I'm guessing Mod(b,n)^d comes to 0 but I'm unsure.