MATLAB: Script doesn’t seem to be execute properly

functioniterationloopsSymbolic Math Toolboxtimewhile loop

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.

Best Answer

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.

Related Solutions

MATLAB: Function to find the next prime number…

Method 1 This functionality does not run in MATLAB and can be used only via MuPAD Notebook Interface.

To create an MuPAD interface use the following code

 mupad
nb = allMuPADNotebooks

Then a interface screen shall pop up where you can proceed by using the nextprime(num) function.

Method 2

You could also create a user defined function to compute the next prime number. This function takes only a non-negative integers as an argument

 function p = nextprime(n)
if (isprime(n))
    p=n;
else
    while(~isprime(n))
        n=n+1;
    end
    p=n;
end
end

MATLAB: Isprime function seems to have poor performance

First, using isprime to test for primality of a Mersenne number is a bad idea. The Lucas-Lehmer test, is as I recall, much more efficient here, and it gives a statement of primality, NOT a statement of isprobable primality. The Miller-Rabin test is an isprobable test, not a proof of primality! That is why multiple internal runs are employed. Do some reading about how these methods can fail.

But more imprtantly, read about Lucas-Lehmer.

https://en.wikipedia.org/wiki/Lucas–Lehmer_primality_test

In fact, sym/isprime has greatly improved over previous ersions. I recall testing it some years ago, and it was pretty slow. But that is not so true now. This means if you have an old release, it might gain you to upgrade. By the way, if I do a test for primality of a number with on the order of 70000 digits, it is a matter typically of 20 minutes, as I recall. (Its been a month or so since I was running massive overnight jobs, searching for prime members of a specific family of primes, thus quasi-modified Woodal primes of a specific class. They can be exceedingly rare if you choose the right family.)

Anyway, you don't want to use isprime to test for a Mersenne prime. Use Lucas-Lehmer. It is REALLY fast in comparison. For a quick implementation...

function s = Lucas_Lehmer(p)
  % Returns true for prime 2^p-1
  
  % don't even bother if p is not prime
  if ~isprime(p)
    s = false;
    return
  end
  
  s = 4;
  M = sym(2)^p - 1;
  for i =  1:p-2
    s = mod(s*s - 2,M);
  end
  s = logical(s == 0);

Seriously, you can't do something much simpler.

Lucas_Lehmer(607)
ans =
  logical
   1
   
isprime(sym(2)^607-1)
ans =
  logical
   1
   
timeit(@() Lucas_Lehmer(607))
ans =
           0.4162580338225
tic,isprime(sym(2)^607-1),toc
ans =
  logical
   1
Elapsed time is 9.295827 seconds.
Lucas_Lehmer(613)
ans =
  logical
   0
% A quick check of the list of known Mersenne primes tells me that 21701 is...
tic,Lucas_Lehmer(21701),toc
ans =
  logical
   1
Elapsed time is 17.560138 seconds.

2^21701-1 is not really that large of a number, but I expect that isprime will take considerably longer. (My own test right now is still running after about 10 minutes.)

tic,isprime(sym(2)^21701-1),toc
???

Finally, if your goal is to find seriously large Mersenne primes, MATLAB is probably not the tool to use if you are looking for primes with millions of digits at some point. But I have found it to work reasonably well in my own play time, searching for primes with 50000-100000 decimal digits. For any seriously large investigation, I'd suggest looking at GIMPS, which seems to be the state of the art.

If your goal is to search for primes in some other family, then there are sometimes tricks to make things more efficient, but they are often strongly dependent on the number family itself. For example, Cullen and Woodall number families have some neat things you can do, to avoid calling isprime too often.

Best Answer

Related Solutions

MATLAB: Function to find the next prime number…

MATLAB: Isprime function seems to have poor performance

Related Question