# Script doesn't seem to be execute properly

3 views (last 30 days)
Manuel Barros on 7 Dec 2018
Commented: Manuel Barros on 10 Dec 2018
I was hoping to acquire some help on how to make my program work efficiently and not take a substantial amount of time to finish:
clear variables
a=...;
p=nextprime(a);
count=0;
limit=200000;
tic
while isprime((p-1)/2)~=1
a=a+1;
p=nextprime(a);
count=count + 1;
if count>limit
break
end
end
toc
This program outputs a number p greater than a such that p is prime and (p-1)/2 is prime. However I've noticed that for any number a greater than approximately 15 digits, the program will take an absurd amount of time to finish, which isn't ideal since I need to test numbers of the order 10^50.

Walter Roberson on 7 Dec 2018
Beyond about 4E15 the distance between adjacent representable doubles becomes greater than 1. p becomes forced to be even (and so not a prime) and p-1 becomes the same as p .
You can do marginally better by switching to uint64, which gets you to about 1.8E19 . But you cannot get beyond that using ordinary numeric forms.
You need to switch to a variable precision toolbox, such as Symbolic Toolbox, or John D'Errico's File Exchange contribution for variable precision integers.

Show 1 older comment
Walter Roberson on 8 Dec 2018
prime tests for a series of numbers are sometimes more efficient as you test more values at one time . It depends on the implementation .
Is the question to find the first such number greater than a or to find some such number greater than a?
Manuel Barros on 8 Dec 2018
It is to find the first prime greater than a such that (p-1)/2 is also prime
Walter Roberson on 8 Dec 2018
Then it is going to depend upon the quality of implementation of isprime() or nextprime() . There is a possibility that it might be faster to test
test_vals = p : 2 : p + 10000;
But that is going to depend on how the isprime() and nextprime() are implemented in the symbolic package.

Christopher Creutzig on 10 Dec 2018
Edited: Christopher Creutzig on 10 Dec 2018
In your code, you spend a lot of time computing the same prime over and over again. Do not start the search at a+1 for the second search, but start after the prime you already found.
It might also be marginally faster to look for the next prime q starting at a/2 such that p=2*q+1 is also prime.
>> tic
>> a = sym('12345678901234567890');
>> q = nextprime(fix(a/2));
>> while ~isprime(2*q+1), q = nextprime(q+1); end
>> toc
Elapsed time is 4.304840 seconds.
>> [q, 2*q+1]
ans =
[ 6172839450617290091, 12345678901234580183]

Stephen Cobeldick on 10 Dec 2018
+1 impressively fast for symbolic math!
Manuel Barros on 10 Dec 2018
Yes thank you, I happened to notice this underlying issue a while after. :)