觀察以下素數表，不能被小於素數n的素數整除的最小數為下一個最小素數。換句話說所謂素數可以看做是座標中不斷使用更小素數的倍數填充後剩餘的最小數。

例：

1】2是最小素數，使用2的倍數不斷填充座標後，未被填充的最小數就是3，於是3就是大於2的下一個素數

2】然後使用3的倍數不斷填充座標，剩餘的未被填充的最小數就是5（4已經被2的倍數填充），可以發現5是大於3的下一個素數

3】依次類推可以順序的求解所有素數

...

隱隱約約覺得，通過以上規律可以推出某範圍內素數個數與素數值的關係，我自己就沒有再推了，有興趣的可以嘗試一下，如果有結果瞭望能分享給我一下，感謝。

然後“talk is cheap，give me the code（廢話少說，放碼過來）”

clc
clear all
close all

numEnter = 100000;
numVector = zeros(1,numEnter);
numVector(1,1) = 1;
disp('2');
divNum = 2;
tempVector = find(~numVector);

tic
while(divNum <= numEnter)
    for iloop1 = 1:length(tempVector)
        % 小於divNum的數字不需要重複判斷，相較於從2開始計算，可以節省1/3的時間（numEnter = 100000時）
        if(mod(tempVector(iloop1),divNum) == 0)
            numVector(1,tempVector(iloop1)) = 1;
        end
    end
    
    tempVector = find(~numVector);       % 尋找所有原始向量中的零元素座標
    if (min(size(tempVector))>0)
        divNum = tempVector(1,1);           % 剩餘的零元素的第一個座標就是下一個素數
    else
        disp('end');
        break;
    end
    disp(divNum);
    % disp 命令會增加耗時5%左右
end
toc
 

以上程式碼的最大特點就是使用排除法對整個數列進行篩選，一旦某個數字被確定是非素數，該數字就會被排除在下一輪的判斷之外，同時一旦被確定了這個數是素數，那麼所有不大於這個數的數字都會被排除在下一輪的判斷之外，執行耗時如下，這是當前已知在不對任何一個素數漏判的前提條件下所需計算量最小的演算法。

Elapsed time is 3.513769 seconds.

掃描驗證一下n = 10000以內範圍中，素數個數與n/lnn之間的關係（掃描時候用倍頻程，應該得到的曲線會更漂亮一些，而且需要的計算量會更少一些）。

clc
clear
close all

range = 1000;
range1 = zeros(1,range);
nlnnVector = zeros(1,range);
primeNumberAcountVector = zeros(1,range);
radioVector = zeros(1,range);

hbar = waitbar(0,'計算進度');
for iloop1 = 1:range
    range1(iloop1) = iloop1*10;
    nlnnVector(iloop1) = range1(iloop1)/log(range1(iloop1));
    
    primeNumberAcount = primeNumberSolute(range1(iloop1));
    primeNumberAcountVector(iloop1) = primeNumberAcount;
    
    radioVector(iloop1) = primeNumberAcountVector(iloop1)/nlnnVector(iloop1);
    
    waitbar(iloop1/range);
end
close(hbar)

figure
subplot(2,1,1)
hold on
plot(range1,nlnnVector)
plot(range1,primeNumberAcountVector);
legend('nlnnVector','primeNumberAcountVector')
hold off
subplot(2,1,2)
semilogx(range1,radioVector)
 

結果如下圖

圖2是兩個結果之間的比值。