% Jake Bobowski
% August 23, 2017
% Created using MATLAB R2014a

% This script calls a function, 'primes_adv.m', to find all of the prime numbers between 'a' 
% and 'b'.  It's an "advanced" version of the function 'primes.m'.  There's
% nothing different about how we're calling the function.  IT's just that
% what's executed within the actual function is slightly more sophisticated
% than what's in 'primes.m'.

% 'primes_adv.m' maintains a file of the sequential prime numbers that it
% has previously found.  Then, if a user asks for prime numbers that it has
% previously found, it just consults the list and tells the user the
% numbers they want to know without doing any searching.  Only if the user
% has asked for prime numbers that it has never found before, does it
% actually implement a search.

% In addition, 'primes_adv.m' uses this list of prime numbers to
% intelligently search for new prime numbers. To check if n is a prime
% number, it is only necessary to check if it is a factor of any of the prime
% numbers less than n.  For example, we don't need to check if 27 is a
% factor of n.  If it is, then 3 is also a factor of n and we would have
% already determined that n is not prime.  By eliminating as many factors
% from our search as possible, we can greatly enhance the efficiency of our
% search for new prime numbers.

% As of Aug. 22, 2017, all of the prime number between 1 and 1e8 have been 
% found and sequentially listed by 'primes_adv.m'.  There are 5,761,455 
% prime numbers less than 100,000,000.

% Note that this is no where close to the largest prime number ever found.
% As I type this, the largest know prime number is 2^(57,885,161) - 1.  A
% number that is a ridiculous 17,425,170 digits long.  However, people do
% not know the sequential list of primes up to this number.

% Note that, after creating this script, I read that people use much faster
% sieve programs to find sequential prime numbers. 

clearvars
format long

% The function that this script calls is named 'primes_adv.m'. It finds all 
% of the prime numbers between input integers a and b.  The code saved in 
% 'primes_adv.m' is shown below without any explanations.  If you open the 
% file 'primes_adv.m', you will find some additional comments that help 
% explain the purpose of the various lines of code.

% To run this function, you will need to have the files 'primes_adv.m',
% 'primes_max.txt', and 'primes_list.txt' in the same directory as the
% script that calls 'primes_adv.m'.  'primes_list.txt' is inside 
% 'primes_list.zip' which can be downloaded from my website.  With a list
% all prime numbers less than 1e8, the zip file is about 14 Mb.  When 
% unzipped, 'primes_list.txt' is about 50 Mb.


% function f = primes_adv(a, b) 
%     fclose('all');
%     if rem(a, 1) == 0 && rem(b, 1) == 0 && a > 0 && b > 0 && b > a 
%         primeMax = dlmread('primes_max.txt');
%         primeNums = dlmread('primes_list.txt');
%         if b <= primeMax
%             i = 1;
%             z = 1;
%             while i <= length(primeNums) && primeNums(i) <= b
%                 if primeNums(i) >= a
%                     piN(z) = primeNums(i);
%                     z = z + 1;
%                 end
%                 i = i + 1;
%             end
%             f = piN;
%         elseif b > primeMax
%             counter = 0;
%             fileID = fopen('primes_list.txt','a');
%             fmt = '%d\n';
%             for i = (primeMax + 1):b
%                 cnt = 0; 
%                 if (i ~= 1 && mod(i, 2) == 1 && mod(sum(int2str(i)-48), 3) ~= 0) && mod(i, 10) ~= 5 || i == 2 || i == 3 || i == 5 
%                     j = 2;
%                     while primeNums(j) < i/primeNums(j - 1) && cnt == 0 
%                         if mod(i, primeNums(j)) == 0 
%                             cnt = 1;
%                             primeMax = i;
%                         end
%                         j = j + 1; 
%                     end
%                     if cnt == 0 
%                         primeMax = i;
%                         primeNums(length(primeNums) + 1) = i; 
%                         dlmwrite('primes_max.txt', primeMax, 'precision', 12, 'delimiter', '\t')
%                         % fileID = fopen('primes_list.txt','a');
%                         % fmt = '%d\n';
%                         fprintf(fileID, fmt, i);
%                         % fclose(fileID);
%                     end
%                 else
%                     primeMax = i;
%                     % dlmwrite('primes_max.txt', primeMax, 'precision', 12, 'delimiter', '\t')
%                 end
%             end
%             fclose(fileID);
%             dlmwrite('primes_max.txt', primeMax, 'precision', 12, 'delimiter', '\t')
%             i = 1;
%             z = 1;
%             while i <= length(primeNums) && primeNums(i) <= b
%                 if primeNums(i) >= a
%                     piN(z) = primeNums(i);
%                     z = z + 1;
%                 end
%                 i = i + 1;
%             end
%             f = piN;
%         end
%     else
%         f = 'ERROR: a and b in primes(a, b) must be positive integers with b > a.'; 
%     end
% end


datestr(clock)
n_start = 1316123820;
n_max = 1316124000;
p_vect = primes_adv(n_start, n_max);
datestr(clock)
length(p_vect)
% A list of the first 5,761,455 sequential prime numbers in ~25 seconds.
