|
何冬州的百度空間Blog 注： wsktuuytyh 即何冬州 三字的五筆編碼
素數(shù)計數(shù)公式全面拉丁化改寫-小有改進-Meissel公式-梅塞爾-Lehmer公式-萊梅=勒梅爾-篩法三種形式-孟慶余公式 2012年06月13日 星期三 6:17 本文標題： 素數(shù)計數(shù)公式全面拉丁化改寫-小有改進-Meissel公式-梅塞爾-Lehmer公式-萊梅=勒梅爾-篩法三種形式-孟慶余公式 前言: 約在2002~2003年我在家閑呆了段時間,專心自學數(shù)論,矩陣論等資料,自己有些心得。 其中有一件小事,就是計算萬以內(nèi)的素數(shù)的個數(shù). 當時自行推導幾個公式,近日發(fā)現(xiàn),與Legendre公式, Meissel公式一致,但沒有Meissel公式這樣簡明.近日看到Meissel公式,與自己的想法一結(jié)合,發(fā)現(xiàn)可以簡化Meissel公式. 而 Lehmer公式,我當時有這樣的思路,作了這樣的計算,不過沒有整理成公式. 此外,看到孟慶余先生的公式,覺得很奇特,竟然使用四舍五入,我還沒有驗證. 而王曉明素數(shù)計數(shù)公式,我在2003年的數(shù)學筆記中已經(jīng)有這樣的計算,但是,我并不覺得有什么特別,因為黎曼和歐拉等人比我走得老遠.我當前試圖精確控制誤差,找到這個公式的小數(shù)部分的誤差限,結(jié)果沒有成功.網(wǎng)上盛傳的一些數(shù)論難題的所謂證明, 所謂公式,都是在拿誤差限不明確的公式在論證,我覺得十分有失嚴謹.當然, 如果公式的誤差限找到,這些證明的過程是可以參考利用的,意義是有的,但是,現(xiàn)在無法肯定. 因此,我提醒數(shù)論愛好者們,包括王曉明先生,孟慶余先生,彎國強先生,多多研讀經(jīng)典,治學嚴謹一些,謙虛務(wù)實,不可累于虛名.雖然我這些話并不是什么新內(nèi)容新道理,但是我強調(diào)一下,亦作為自勉吧. 而王曉明先生的索數(shù)普遍公式,實際是構(gòu)造了質(zhì)數(shù)連乘積的縮系(既約剩余系,簡化剩余系),我于2003年受臺爾曼公式和閩嗣鶴·嚴士健《初等數(shù)論》的啟發(fā), 給出了構(gòu)造縮系的更簡明的方案,例如 ( ( 1/2+{1,2}/3+{1,2,3,4}/5+{1,2,…,6}/7 ) mod 1 ) *2*3*5*7 其中a mod 1,即是a的小數(shù)部分. . 2004-2005年,受洪伯陽先生的《數(shù)學寶山上的明珠》的啟發(fā),在研究中國剩余定理有所改進之后,給出了另外的表示.2010年,看過了王曉明先生的文章后,發(fā)現(xiàn)始終沒有解決誤差限問題和算法復雜性問題,所以,我曾經(jīng)致信王曉明先生, 認為素數(shù)普遍公式,并沒有完成,也不成熟. 因此,我略整理了一下網(wǎng)上搜集來的經(jīng)典結(jié)論,打算從中獲得一些啟發(fā),鞏固自己的數(shù)論基礎(chǔ). 其中有我個人的心得,如有不足不妥,請指教.謝謝. AAAMeissel公式及其小小改進 以下說Meissel公式之我的改寫,結(jié)合了我舊日的一些思路,有空我會查一下舊日的資料,或者能發(fā)現(xiàn)一些或啟發(fā)一些好結(jié)果.歷史總會發(fā)展的.古人也是能超越的,包括過去的自己. Y(m,n)=Phi(m,n)=card{i; i<=m, gcd(i, 2*...*p(n))=1)注: card表示集合元素個數(shù)(集合的基數(shù)或勢) 即m以內(nèi)與2,…,p(n)互質(zhì)的數(shù)的個數(shù),與歐拉函數(shù)相似,故采用Phi記號,希臘字母是Φ.這里p(n)則是prime(n),第n個素數(shù).在我的數(shù)論筆記中,我稱它為泛縮系計量. Y(m,n)=Y(m,n-1)-Y([m/p(i)],n-1)注：這是個很簡單的遞推公式。注意,顯然,
Y(m,0)=m,因為所有數(shù)均與1互質(zhì). 此外還有下面的性質(zhì): Y(a*prod(2,…,P(k))+r, k)=a*Y(prod(2,...,p(k))+Y(r,k)=a*(2-1)(3-1)...(p(k)-1)+Y(r,k) 其中prod表示連乘積,即Π(p(i)),素數(shù)p(i)的序號由1取到k. 注: 我把這個式子也用于計算自變量為負數(shù)的情況,例如,定義Y(a*prod-r, k)= a*Y(prod(2,...,p(k))+Y(-r,k),計算也此一致,簡化方案是利用 [-x]=-ceil(x),亦作ceiling(x),例如ceil(2.3)=3, ceil(2)=2,向上取整. 記c=pi(m^(1/3)), s=pi(m^(1/2)),這里 pi(x)表示x以內(nèi)的素數(shù)個數(shù),即π(x).注意,pi(c)<=cbrt(m), pi(s)<=sqrt(m). 即c=pi(cbrt(m)),cuberoot(立方根)以下的素數(shù)個數(shù). s=pi(sqrt(m)), squareboot(平方根)以下的素數(shù)個數(shù).易記. 則pi(m)=Y(m,c)+(s(s-1)-(c-1)(c-2))/2-sum{pi(m/p(i))}{sum表示求和,即Σ.這里范圍為i=c+1到s} 注:這個公式比Meissel公式形式簡潔,計算容易了一點點,應(yīng)當算是略有改進吧。 Φ(m,n)=Φ(m,n-1)-Φ([m/p(i)],n-1) 給定一個自然數(shù)m,若n=π(m^(1/3))且μ=π(m^(1/2))-n,那么： π(m)=Φ(m,n)+n(μ+1)+(μ^2-μ)/2-1-Σπ(m/p(n+k)) ,求和限:k={1,…,μ}.
例:計算m=1000時的pi(m). 易算得Y(1000, 4)=Y(1000,3)-Y(142,3)=266-38 =228 注：Y(1000,3)的計算，參見我的一個答題(題中用的是類似Legendre的算法)，見CCC相關(guān)題 s=card {2, 3, 5 ,7 ,11, 13, 17 ,19, 23 ,29, 31 } =11, c=pi (10) =4 這里card表示集合元素數(shù)量(基數(shù)或勢). 1/2*(s(s-1)-(c-1)(c-2)) =55-3=52 序列 { pi (m / { 11, 13, 17, 19, 23, 29, 31} ) } = { pi (90, 76, 58, 52, 43, 34, 32) } 注意:首項為 25 -4,反差分(相鄰項前項減后項)序列為 {3, 5, 1, 1, 3, 0},即相鄰兩數(shù)間的素數(shù)個數(shù) ={ 24, 21, 16, 15, 14, 11, 11 } 其和為 112 故 pi(1000) =228 + 52 -112 =168.
AAAAAA篩法淺說 1一般篩法,先用少數(shù)幾個素數(shù)篩去合數(shù),即與它們的乘積不互素的數(shù).然后再用更大的素數(shù),在剩下的數(shù)里面類似篩去.如厄拉多塞篩法. 2有些篩法,在剩下的數(shù)中篩選時,根本不用到更大的素數(shù),只是原來的素數(shù)的縮系中的同余類 (我稱為泛縮系)來構(gòu)造出剩下來的數(shù)中的所有合數(shù),排除它們就得到質(zhì)數(shù).這種篩法, 可稱為反篩法,最初使用的是辛答拉姆.其實原理很簡單,就是從奇數(shù)中篩去奇合數(shù),而奇合數(shù)由 集合 {2n+1; n>=1}元素的乘積構(gòu)成. 當初辛答拉姆新開一路,的確可貴.但是,辛答拉姆當初將奇合數(shù)全部 折半取整, 在判別素數(shù)時又乘2加1,其實是在繞彎子,并不是必要手段,也沒有任何神秘性可言,好處是可以減小存儲或說明該數(shù)表的構(gòu)造可以不依賴于乘積. 個人認為,不必過于在這個彎上繞來繞去,要看破這一點. 3就是下面要說到的Lehmer篩法.當然得算是一種新的篩法,思路已經(jīng)又有所不同了,不應(yīng)當算是厄拉多塞篩法的改進.不然, 那個篩法不是源于厄拉多塞篩法這個歷史源頭呢?因此,我專稱為Lehmer篩法. 4.我曾經(jīng)考慮用素數(shù)的平方為篩子,篩得所有的無平方因子數(shù),然后再在其中篩去合數(shù).我發(fā)現(xiàn),無平方因子數(shù)的分布規(guī)律也很復雜.但是,我覺得研究起來可能會容易一些.
BBB Lehmer篩法與Lehmer素數(shù)計數(shù)公式. Lehmer的算法,實際是另一種篩法.我也這樣想過,這樣算過.想不到是閉門造車,而且是舊車和小車,竟沒有造出古人造過的大車. Lehmer用他的方法算到了pi(10^10),即100億內(nèi)的素數(shù)的個數(shù).我僅僅是此這種思路驗證過 pi(10000) ,pi(1000), 具體情況,還得找我舊日的筆記. 于是,這里只介紹一下思路,有興趣的朋友,可在網(wǎng)上搜索 素數(shù)計數(shù)函數(shù) OR 篩法計算公式 即可找到相關(guān)資料. Lehmer的篩法,真正是篩去合數(shù),根據(jù)合數(shù)的因子個數(shù)來篩去.定義 P_k(m,n)=card {i; i<=m, i有k個素因子,且素因子>prime(n)} 其中prime(n)也作p(n),指素數(shù)列的第n項. 注:其中這k個素因子可能相等也可不等. 顯然,根據(jù)數(shù)的質(zhì)因子的個數(shù)為0,1,2,3, ... ,可以將所有正整數(shù)進行分類.即 m =P_0+P_1+P_2+...+P_∞, 而對于與2*3*…*p(n)互質(zhì),即只有>p(n)的素因子的數(shù),其個數(shù)為 Y(m,n)=P_0(m,n)+P_1(m,n)+P_2(m,n)+… 這里的Y(m,n)同前面講Meissel公式時的相同,而lehmer給出了另一種算法. 注意,在一定范圍內(nèi)的數(shù),素因子個數(shù)是有限的,即上式的后面從某個項開始,全部為0.例如, Lehmer就考慮 m的立方根以下的最大素數(shù)p(c), m以內(nèi)有三個以上比這個素數(shù)還大的素因子的數(shù),顯然是沒有的.因而就不必再往上計算了.即下面所講到的,n>=c時,P_k(m,n)=0, k>=3 下面開始計算. P_0(m,n)=1, 注:顯然可以直觀地看出.即沒有素因子的數(shù),當然只有 {1}一個. P_1(m,n)=pi(m)-n 注:有一個大于p(n)的素因子的個數(shù),即大于p(n)的素數(shù)的個數(shù). P_k(m,n). k>=3, Lehmer設(shè)法使其值為0,從而不必計算.條件是n>=c=pi(m^(1/3)).即p(n)>=p(c).注意,是>=,大于等于.一個相關(guān)細節(jié):cbrt(m)>=p(c), p(c)是<=cbrt(m)=m^(1/3)的最大素數(shù). P_2(m,n)的計算.顯然,當n>=s=pi(m^(1/2)),即p(n)>=sqrt(m)=m^(1/2)時, P_2(m,n)=0.注意,是>=,大于等于.因此,計算P_2(m,n)時,取n<s,只有與此對應(yīng)的質(zhì)數(shù)參加主要計算過程. 綜上,當c<=n時,Y(m,n)=1+pi(m)-n+P_2(m,n). 且P_2(m,n)的計算主要與c<=n<s對應(yīng)的質(zhì)數(shù)相關(guān).即p(c)<=p(n)<p(s).一個相關(guān)細節(jié): p(s)是<=sqrt(m)的最大素數(shù). Lehmer公式的重點: P_2(m,n)=sum {pi(m/p) -pi(p)+1} {求和范圍: p(n)<p<=p(s) } 何冬州注: =sum{pi(m/p)}-sum{n,n+1,…,s-1}=sum{pi(m/p)}-(n+s-1)(s-n)/2 這里的n值可在c<=n內(nèi)選取.為了計算方便,還可限定在c<=n<s內(nèi)選取,這一點較Meissel公式靈活,因此有些資料上也說是Meissel公式的推廣.其實不僅僅是推廣,而是一種新篩法思路.下面會舉例,取n=5,也可計算pi(1000). 于是, pi(m)=Y(m,n)-1+n-P_2(m,n), c=<n<s. 何冬州注:即是pi(m)=Y(m,n)-sum{pi(m/p)}+n-1+(n+s-1)(s-n)/2 我的簡化之一為:pi(m)=Y(m,n)-sum{pi(m/p)}+(s(s-1)-(n-1)(n-2))/2 僅僅是使其易記易用了一點點.但對于我個人而言,很有好處.當取n=c時,與Meissel公式一致. { 其中Y(m,n)的計算過程,可以用前面講 Meissel公式時用到的,這里重錄一下: Y(m,n)=Y(m,n-1)-Y([m/p(i)],n-1)注：這是個很簡單的遞推公式。注意,顯然, Y(m,0)=m,因為所有數(shù)均與1互質(zhì). 此外還有下面的性質(zhì): Y(a*prod(2,…,P(k))+r, k)=a*Y(prod(2,...,p(k))+Y(r,k)=a*(2-1)(3-1)...(p(k)-1)+Y(r,k) } Lehmer公式應(yīng)用舉例: 用來計算pi(m),m=1000, c=4, s=pi {2,3,5,7,11,13,17,19,23,29,31}=11 注:其中利用到了前面算過的結(jié)果Y(1000,4)=228.另外,Y(90,4)=Y(90,3)-Y(12,3)= card {13到90以內(nèi)形如30k+{1,7,11,13,17,19,23,29}的數(shù)} =21.(=5+16或24-3) sum{ pi(m\{13,…,31})}=112-pi(1000\11)=112-pi(90)=112-24=88.注:這里利用了前面的計算結(jié)果sum{ pi(m\{11,…,31})}=112. 11*10/2-4*3/2=55-6=49 于是pi(1000)=207+49-88=168 我們可以看到,Lehmer公式可以充分轉(zhuǎn)化與利用已有成果,我想lehmer必有一套流程,可以在特定的一組數(shù)內(nèi)批量利用,計算pi(x)值. 我曾經(jīng)在不知Meissel公式與Lehmer公式時獨立計算pi(x),也曾有此想法,作了大量筆記.
CCCAAA：Lengdre素數(shù)計數(shù)公式,容斥原理=包含排除原理=逐步淘汰原則 參考: 歐拉函數(shù),既約剩余系=簡化剩余系=縮系,縮系的構(gòu)造, Mobius函數(shù) 基礎(chǔ)：記號同前： Y(m,n)=Phi(m,n)=card{i; i<=m, gcd(i, 2*...*p(n))=1)注: card表示集合元素個數(shù)(基數(shù)或勢) 即m以內(nèi)與2,...,p(n)互質(zhì)的數(shù)的個數(shù),與歐拉函數(shù)相似,故采用Phi記號,希臘字母是Φ. 記s=pi(m^(1/2)),這里 pi(x)表示x以內(nèi)的素數(shù)個數(shù),即π(x).注意, pi(s)<=sqrt(m);而s=pi(sqrt(m)),即 squareboot(平方根)以下的素數(shù)個數(shù).易記. 易見 Y(m, s) = pi (m) -pi (s) +1,即得 Lengdre素數(shù)計數(shù)公式: pi(m) = -1 + pi(s) + Y(m, s)
以m=100為例說明,為方便,用a \ b 表示 a/b的整數(shù)部分.此時pi(s)=pi(sqrt(100))=pi(10)=4 Y(m,s)=card{i; i<=100, gcd(i, 2*...*7)=1) = 100 -1- card {100內(nèi)能被2, 3, 5, 7之一整除的數(shù)} =100 - 100 \ { 2, 3, 5, 7 } + (100 \2) \ {3, 5, 7} + (100\3) \ {5,7} + (100 \5) \ 7 -(100 \2\3 ) \ {5,7} -(100\2\5) \7 =(100 - 50-33-20-14 + 16+10+7 + 6+4 +2 -3-2 -1)=22 故pi(100)=-1+4+22=25
注:可以看到,有些項是可以繼承利用的;列成表格,類似于卡諾圖,元素個數(shù)構(gòu)成2的方冪. 注:試利用Y(m,n)=Y(m,n-1)-Y([m/p(i)],n-1)簡化計算: 100-14-(100 \{2,3,5}-14\{})+(100\{2*3,2*5,3*5}-14\{})-(100\(2*3*5)-14\())=(100-14-(50+33+20-7-4-2)+(16+10+6-2-1)-(3))=22 注:試利用 Meissel公式來計算,看看改進了什么.快在什么地方. pi(m)=Y(m,c)+(s(s-s)-(c-1)(c-2))/2-sum{pi(m/p(i))} {sum表示求和,即Σ.這里范圍為i=c+1到s} m=100, c=card {2,3}=2, s=4, pi(100)=Y(100,2)+ (4*3-1*0)/2-sum{pi {100/5, 100/7}} =(100-50-33+16)+6-sum(card{2,3,5,7,11,13,17,19}, card{2,3,5,7,11,13})=17+22-8-6=25.簡單多了.
CCC相關(guān)題: 1到1000中所有不能被2、3、5整除的自然數(shù)有多少個？算式怎么列？ http://zhidao.baidu.com/question/435220655.htm 其中的容斥原理,就是Legendre算法的原理. 解：以下用[]表示取整數(shù)部分,即高斯取整函數(shù). 計算過程是將下面{}內(nèi)的各行的計算項取代數(shù)和. { 1000, -[1000/2]-[1000/3]-[1000/5] = -500 -333-200 =-1033 +[[1000/2]/3] +[[1000/2]/5] +[[1000/3]/5] =[500/3]+[500/5]+[333/5]=166+100+66=332 -[ [[1000/2]/3]/5 ]=-[166/5] =-33 } =1000-(500+333+200) +(166+100+66)-33 =-33+332-33=266 原理:容斥原理 (又稱作包含排除原理或逐步淘汰原則,可參考百度百科-容斥原理) 推廣: 1到X中所有不能被2, 3, 5, ... , p(i), ..., p(n)整除的自然數(shù)的個數(shù) =sum{ +X, -{ [X/2]+[X/3]+...+...+[X/p(n)] } +{ +[X/2/3]+[X/2/5+....+[X/2/p(n)] + [x/3/5]+... +[x/2/p(n)] + ... +... + +[X/p(n-1)/p(n)] } -... +... + (-1)^n * [X/2/3/.../p(n)] } 注#1:這里的[]表示取整數(shù)部分,即高斯取整函數(shù). 注#2:計算時可以使用 [X/(2*3)]=[ [X/2] /3 ]= [ [X/3] /2 ] 注#3:在數(shù)論中,歐拉函數(shù)有類似的形式與性質(zhì) 注#4:代數(shù)和的正負號=(-1)^數(shù)的奇次方質(zhì)因子的個數(shù)(或者說,去除平方因數(shù)之后的剩下質(zhì)因子的個數(shù)為奇則取負,為偶數(shù)包括0則取正).這在數(shù)論中稱作mobius函數(shù)(莫比烏斯),參考百度百科-數(shù)論函數(shù)-積性函數(shù). |
|
|
