費(fèi)馬數(shù)乘積形式與平方和形式轉(zhuǎn)換研究

周利榮 胡天磊

摘要：該文分析了費(fèi)馬數(shù)素因子乘積形式轉(zhuǎn)換為平方和形式的算法，費(fèi)馬數(shù)平方和形式轉(zhuǎn)換為素因子乘積形式的算法。利用轉(zhuǎn)換算法2在極短的時(shí)間內(nèi)（0.02 s）將62位素?cái)?shù)934616397153579777691635581996068965840512375416381885802 80321分解成平方和形式。綜合利用兩種算法部分分解費(fèi)馬數(shù)F12。

關(guān)鍵詞：費(fèi)馬數(shù);費(fèi)馬平方和定理;呂卡定理

中圖分類號(hào)：TP312? ? ? 文獻(xiàn)標(biāo)識(shí)碼：A

文章編號(hào)：1009-3044（2021）28-0114-03

開放科學(xué)（資源服務(wù)）標(biāo)識(shí)碼（OSID）：

1 背景

表達(dá)式Fn=[22n]+1稱為費(fèi)馬數(shù)，當(dāng)n取0、l、2、3、4時(shí)的值分別為3、5、17、257、65537，并均為素?cái)?shù)。由此，費(fèi)馬提出一個(gè)猜想：形如Fn=[22n]+l的數(shù)一定為素?cái)?shù)。數(shù)學(xué)家歐拉發(fā)現(xiàn)，當(dāng)n取5時(shí)，F(xiàn)5=4294967297=641*6700417不為素?cái)?shù)。費(fèi)馬素?cái)?shù)除了被費(fèi)馬本人所證實(shí)的那五個(gè)外竟然沒有再發(fā)現(xiàn)一個(gè)。費(fèi)馬數(shù)問題促進(jìn)素性判定算法和整數(shù)的分解算法的研究。

2 費(fèi)馬數(shù)性質(zhì)

有關(guān)定理及性質(zhì)：

1）費(fèi)馬平方和定理：任何一個(gè)4p+1型素?cái)?shù)可表成兩個(gè)整數(shù)的平方和，且表示方法是唯一的。

2）呂卡定理：費(fèi)馬數(shù)Fn=[22n]+1的素因數(shù)必具p=k.[2n+2]+l的形式，其中k為正整數(shù)[1]。

3）性質(zhì)1：當(dāng)n≥2時(shí)，費(fèi)馬數(shù)Fn=[22n]+1的末位數(shù)為7[2]。

4）性質(zhì)2：（[u2+v2]）*（[a2+b2]）=[（ua+vb）2] +（[va-ub）2]。

5）性質(zhì)3：（[u2+v2]）*（[a2+b2]）=[（ua-vb）2] +（[va+ub）2]。

3 費(fèi)馬數(shù)素因子乘積形式轉(zhuǎn)換成平方和形式

3.1 轉(zhuǎn)換的依據(jù)

由呂卡定理可知費(fèi)馬數(shù)的素因數(shù)必具k.[2n+2]+l的形式，即為4p+1型的素?cái)?shù)。由費(fèi)馬平方和定理可知任何一個(gè)4p+1型素?cái)?shù)可表成兩個(gè)整數(shù)的平方和，且表示方法是唯一的。由性質(zhì)2：（[u2+v2]）*（[a2+b2]）= [（ua+vb）2] +（[va-ub）2]，可知整數(shù)的平方和的乘積仍為平方和。因此，費(fèi)馬數(shù)一定可以轉(zhuǎn)換成整數(shù)的平方和形式。

3.2 轉(zhuǎn)換算法1（窮舉法）

輸入：費(fèi)馬數(shù)的素因數(shù)pi（1≤i≤k）。

輸出：費(fèi)馬數(shù)的平方和形式。

1）i=1;

2）tmp=pi ，t=[pi]-1;

3）tm=tmp-[t2];

4）如果tm 是完全平方數(shù)，轉(zhuǎn)（5），否則 t--，轉(zhuǎn)（3）;

5）保存pi的平方和形式;

6）i++，如果 i≤k，轉(zhuǎn)（2），否則轉(zhuǎn)（7）;

7）由性質(zhì)2輸出費(fèi)馬數(shù)的平方和形式，程序結(jié)束。

設(shè)（x，y）是滿足pi=x2 + y2的格點(diǎn)，則（y，x）也是滿足pi=x2 + y2的格點(diǎn)，即格點(diǎn)在第一象限的分布滿足對(duì)稱性，因此，t的搜索區(qū)間為[[pi/2]，[pi]-1]。

3.3 轉(zhuǎn)換算法2（利用素因子平方和之間的比例關(guān)系）

性質(zhì)2、性質(zhì)3的成立是顯而易見的。由性質(zhì)2、性質(zhì)3可知，兩因子費(fèi)馬數(shù)Fn=p1*p2=（[u2+v2]）*（[a2+b2]）=[（ua+vb）2] +（[va-ub）2]=[（vb-ua）2] +（[va+ub）2]。即兩因子費(fèi)馬數(shù)有兩種平方和形式。由于Fn=[22n]+l=[（22n-1）2+12]，因此，[（22n-1）2+12]必定是其中一種平方和形式。設(shè)u

推論 1：vb-ua =1。

推論 2：[uv]≈[ba]。

推論 3：[uu2+v2]≈[ba2+b2]。

推論 4：[vu2+v2]≈[aa2+b2]。

對(duì)于兩因子費(fèi)馬數(shù)，有了上述推論，如果已將較小因子p1分解成平方和形式p1=（[u2+v2]）。利用推論4可得a的近似值，從a的近似值開始遞增將p2分解成平方和形式[a2+b2]，這樣極大地減少了搜索范圍。

對(duì)于兩因子費(fèi)馬數(shù)Fn=p1*p2， p1

1）用窮舉法將較小因子p1分解成平方和形式p1=（[u2+v2]）;

2）計(jì)算c=v/[p1];

3）由p2=（[a2+b2]）及a/[p2]≈c將較大因子p2分解成平方和形式p2=（[a2+b2]）。

4）由性質(zhì)2輸出費(fèi)馬數(shù)的平方和形式。

例：F5=4294967297=641*6700417。p1=641，p2=6700417。窮舉法分解p1=641=（[42]+[252]），設(shè)u=4，v=25，計(jì)算c=25 /[641 ]=0.987，由推論4知：a/[p2]≈c=0.987得a≈2555，循環(huán)兩次即可分解6700417 =（[25562]+[4092]），設(shè) a=2556，b=409。ua+vb= 4* 2556+25 * 409 =20449，va-ub=25*2556-4*409=62264， F5=（[204492]+[622642]）。

例：F8= 115792089237316195423570985008687907853269984665640564039457584007913129639937=1238926361552897 *93461639715357977769163558199606896584051237541638188580280321。p1=1238926361552897，p2=93461639715357977769163558199606896584051237541638188580280321。窮舉法分解p1=1238926361552897=（[242465592+255153042]）， u=24246559， v=25515304，計(jì)算c=25515304 /[1238926361552897 ]=0.7248998337342965585004679297962，由a/[p2]≈c得a≈7008009763344264886811252498144，循環(huán)兩次即可分解p2 =（[70080097633442648868112524981452]+[66595374368066614409021877834642]）。設(shè)a=7008009763344264886811252498145，b=6659537436806661440902187783464。ua+vb= 2424655 *7008009763344264886811252498145+25515304 * 6659537436806661440902187783464 =33984024439900551177939471112034026611，va-ub=25515304*7008009763344264886811252498145-24246559*6659537436806661440902187783464=17340632172455487023654788790090010704， F8=（[173406321724554870236547887900900107042]+[339840244399005511779394711120340266112]）。

3.4 算法的效率比較

算法2由于利用推論4，算法的時(shí)間主要用于第（1）步：用窮舉法將較小因子p1分解成平方和形式，（2）（3）（4）步所用時(shí)間極少（分解F5的第二個(gè)素因子p2 為平方和形式僅用0.0001s;分解F6的第二個(gè)素因子p2 為平方和形式僅用0.001s;分解F7的第二個(gè)素因子p2 為平方和形式僅用0.01s;分解F8的第二個(gè)素因子p2 （62位素?cái)?shù)）為平方和形式僅用0.02s，是目前被分解的最大素?cái)?shù)）。因此，算法2的效率遠(yuǎn)優(yōu)于算法1。上表中，算法1計(jì)算F7、F8較大因子p2平方和形式所需要時(shí)間=循環(huán)一次時(shí)間*循環(huán)次數(shù)估算得到?？梢姰?dāng)n≥7時(shí)，算法1的效率極低，幾乎不可能。

4 費(fèi)馬數(shù)平方和形式轉(zhuǎn)換成乘積形式（適用于兩因子費(fèi)馬數(shù)）

4.1 轉(zhuǎn)換算法3（解方程組法）

主要利用性質(zhì)2：（[u2+v2]）*（[a2+b2]）= [（ua+vb）2] +（[va-ub）2]。如果已知費(fèi)馬數(shù)的分解形式Fn=x2 + y2，由性質(zhì)2可知，存在u、v、a、b，使得[（ua+vb）2] +（[va-ub）2]= x2 + y2= （[u2+v2]）*（[a2+b2]），只要能求出u、v、a、b，則可將n分解成素因子乘積形式。

要求出u、v、a、b的值，必須有四個(gè)方程組成的方程組。而條件只有兩個(gè)方程：x = ua+vb，y = va-ub。但x，y具有x= x1+x2，y = y2-y1的 形式，因此，可設(shè)x= x1+x2，y = y2-y1，則x2 + y2= （x1+x2）2 + （y2 -y1）2。如果x1，x2，y1，y2滿足x1=[ ua，]x2= vb，y1=[ va]，y2=[ ub]，則計(jì)算u=gcd（x1，y2），v= gcd（x2，y1），則得n的一個(gè)因子（[u2+v2]），分解成功。

如果x1，x2，y1，y2滿足x1=[ ua，]x2= vb，y1=[ va]，y2=[ ub]，則x1x2=uavb，y1y2=vaub，因此y1y2=x1x2，則必有方程組：

y1y2=x1x2……①

y1-y2=y ……②

由于y已知，只要知道x1，x2，則可求y1，y2，進(jìn)而求出u、v、a、b，則可將n分解成素因子乘積形式。由②得：y2=y+y1，代入①得：[y12]+y.y1-x1.x2=0，是一個(gè)一元二次方程，解方程得：y1=（-y+[y2-4x1x2]）/2，y2= y+（-y+[y2-4x1x2]）/2。

推論 5：當(dāng)n≥2時(shí)，費(fèi)馬數(shù)Fn= x2 + y2中x，y必有一個(gè)為奇數(shù)。

由性質(zhì)1不難推出此結(jié)論。設(shè)x為奇數(shù)，如果x1=[x/2];x2=x-x1，則有x2-x1=1。由3.3節(jié)的分析可知：（ua-vb）=1，（[va+ub]）=[22n-1]，又有：ua+vb=x ，va-ub=y 。因此，x1，x2，y1，y2必滿足x1=[ ua，]x2= vb，y1=[ va]，y2=[ ub]，[y2-4x1x2]=[（va-ub）2]+4uavb=[（va+ub）2]是完全平方數(shù)。y1=（-y+[y2-4x1x2]）/2，y2= y+（-y+[y2-4x1x2]）/2必為整數(shù)。因此，算法只要循環(huán)一次即可。

算法3描述如下：

1）求出x1=[x/2];

2）計(jì)算x2=x-x1;

3）解方程組y1y2=x1x2，y1-y2=y;

4）y1、y2必為整數(shù)，由u=gcd（x1，y2），v= gcd（x2，y1）求出u，v，將n分解，程序結(jié)束。

例：F6= 18446744073709551617= x2 + y2=[14387937592+40468032562]，x=1438793759，y= 4046803256，由于x為奇數(shù)，令x1=719396879，x2=719396880， 解方程組y1y2=x1x2，y2-y1=y得y1=（-y+[y2-4x1x2]）/2=124082020， y2= y+y1=4170885276為整數(shù)。u= gcd（y2，x1）=89，v=gcd（y1，x2）=516。（[u2+v2]）=274177，18446744073709551617=274177*67280421310721。

4.2 轉(zhuǎn)換算法4（求最大公約數(shù)法）

由4.1節(jié)的分析可得到如下四個(gè)方程：

vb-ua =1? ? ?……①

va+ub=[22n-1]? ?……②

ua+vb=x? ? ?……③

va-ub=y? ? ?……④

①+③得：2vb= x+1，vb= （x+1）/2;②+④得：2va=[22n-1]+y，va=（[22n-1]+y）/2，v=gcd（（x+1）/2，（[22n-1]+y）/2）。

③-①得：2ua= x-1，ua= （x-1）/2;②-④得：2ub=[22n-1-]y，ub=（[22n-1-]y）/2，u=gcd（（x-1）/2，（[22n-1]-y）/2）。

算法4描述如下：

1）求出v=gcd（（x+1）/2，（[22n-1]+y）/2）;

2）求出u=gcd（（x-1）/2，（[22n-1]-y）/2）;

3）得Fn的一個(gè)因子（[u2+v2]）;

4）將Fn分解n=（[u2+v2]）*[Fn/（[u2+v2]）]，程序結(jié)束。

例：F7=340282366920938463463374607431768211457=

x2 + y2=[163823502215354644792+84794438579364025042]，

設(shè)x=16382350221535464479，y=[8479443857936402504]，

v=gcd（（x+1）/2，（[22n-1]+y）/2）=208648999，u= gcd（（x-1）/2，（[22n-1]-y）/2）=126945596。（[u2+v2]）=59649589127497217，340 2823669209384634633746074317682114577=59649589127497 217 * 5704689200685129054721。

4.3 算法的效率比較

算法3由于第3）步解方程組需要時(shí)間比較多，效率低于算法4，但效率仍然是極高的。從上表可知，對(duì)于兩因子費(fèi)馬數(shù)，在已知費(fèi)馬數(shù)平方和形式的情況下，將其轉(zhuǎn)換成素因子乘積形式的效率是極高的。

5 兩種轉(zhuǎn)換算法在分解多因子費(fèi)馬數(shù)F12中的應(yīng)用

5.1 素因子乘積形式轉(zhuǎn)換成平方和形式

第一步：利用呂卡定理得到費(fèi)馬數(shù)F12的兩個(gè)小素因子。

由呂卡定理可知費(fèi)馬數(shù)的素因數(shù)必具k.[2n+2]+l的形式，k從1開始搜索，當(dāng)k=7時(shí)， p1=7*[214]+1=114689且p1| F12，得到F12的第1個(gè)素因子p1。將F12除以p1，k又從1開始搜索，當(dāng)k=1588時(shí)，p2= 1588*[214]+1=26017793且p2|（F12/ p1），得到F12的第2個(gè)素因子p2?？偤臅r(shí)約為0.03s。

第二步：利用算法1將費(fèi)馬數(shù)F12的兩個(gè)小素因子分解成平方和形式，p1=114689=[2602]+[2172]，p2=26017793=[20472]+[46722]。耗時(shí)約為1.9s。F12是多因子費(fèi)馬數(shù)，推論1-4并不成立，算法2失效。

第三步：利用性質(zhì)2、性質(zhì)3將費(fèi)馬數(shù)F12的兩個(gè)小素因子乘積C13=p1* p2=2983954661377分解成兩種平方和形式，2983954661377=[15460442]+[7705212]，x=770521，y=1546044;2983954661377=[16589192]+[4816042]，n=[481604]，m=1658919。耗時(shí)約為0.02s。

5.2 平方和形式轉(zhuǎn)換成素因子乘積形式

由性質(zhì)2、性質(zhì)3可得到如下四個(gè)方程：

ua-vb =n? ?……①

va+ub=m? ?……②

ua+vb=y? ? ……③

va-ub=x? ? ……④

③-①得：2vb= y-n，vb=（y-n）/2;②+④得：2va=m+x，va=（m+x）/2，v=gcd（（y-n）/2，（m+x）/2）。

①+③得：2ua= y+n，ua=（y+n）/2;②-④得：2ub=m[-]x，ub=（m[-]x）/2，u=gcd（（y+n）/2，（m-x）/2）。

將C13轉(zhuǎn)換成素因子乘積形式，可用如下算法。

算法5描述如下：

1）求出v=gcd（（y-n）/2，（m+x）/2）;

2）求出u=gcd（（y+n）/2，（m-x）/2）;

3）得Fn的一個(gè)因子（[u2+v2]）;

4）將Fn分解n=（[u2+v2]）*[Fn/（[u2+v2]）]，程序結(jié)束。

在已知2983954661377=[15460442]+[7705212]，x=770521，y=1546044;2983954661377=[16589192]+[4816042]，n=[481604]，m=1658919的情況下，v=gcd（（y-n）/2，（m +x）/2）=260，u=gcd（（y+n）/2，（m -x）/2）=217，p1=[u2+v2] =114689。

∴2983954661377=114689*26017793。

即將C13=2983954661377分解成素因子乘積形式，耗時(shí)約為0.02s。

5.3 研究兩種形式轉(zhuǎn)換的意義

費(fèi)馬數(shù)的分解算法有：橢圓曲線分解算法（分解F10、F11）、連分?jǐn)?shù)算法（分解F7）、數(shù)域篩算法（分解F9）、Pollard's rho分解算法（分解F8）[3-4]。由文獻(xiàn)[5]F12已經(jīng)分解出六個(gè)素因子，還有一個(gè)1133位的末尾合數(shù)C1133未分解。即F12 = 114689 * 26017793 * 63766529 * 190274191361 *1256132134125569 * 568630647535356955169033410940867804839360742060818433 * C1133。

上述費(fèi)馬數(shù)的分解算法：橢圓曲線分解算法、連分?jǐn)?shù)算法、數(shù)域篩算法、Pollard's rho分解算法對(duì)C1133顯然是無效的。如果知道C1133的兩種平方和形式，由算法5可求出C1133的一個(gè)素因子p7，且效率極高，利用Miller–Rabin算法[6]、ECPP、APRCL、AKS[7]等算法判定C1133/ p7的素性，如果為C1133/ p7素?cái)?shù)，則F12完全分解，否則求出C1133/ p7的兩種平方和形式，由算法5可求出C1133/ p7的一個(gè)素因子p8，如此循環(huán)直至C1133完全分解。

將C1133這樣的大整數(shù)分解為平方和形式尚無有效算法，如能像算法2一樣找出多因子費(fèi)馬數(shù)的素因子平方和變量（u、v、a、b）之間的比例關(guān)系，由p1、p2、p3、p4、p5、p6求出C1133的平方和形式，是一種可能途徑。這為分解費(fèi)馬數(shù)提供了新的思路和可能性。

C++中用unsigned _int64表示整數(shù)范圍為[0，[264]]，即最大整數(shù)為18446744073709551616，費(fèi)馬F6= 18446744073709551617大于此最大整數(shù)。NTL庫是一個(gè)用于大整數(shù)運(yùn)算的C++庫，可以用于任意長度的整數(shù)計(jì)算及整數(shù)的多項(xiàng)式算法，以上實(shí)驗(yàn)結(jié)果均在win10系統(tǒng)、VS 2017軟件+NTL、intel 酷睿處理器i3? @3.7G hz /8GB內(nèi)存條件下取得。

參考文獻(xiàn)：

[1] 賈耿華.關(guān)于費(fèi)馬數(shù)的研究[D].成都：成都理工大學(xué)，2006.

[2] 劉培杰.費(fèi)馬數(shù)[J].自然雜志，1991，13（8）：608-612.

[3] 周利榮，胡天磊.基于萊梅素?cái)?shù)判定定理的安全素?cái)?shù)構(gòu)造算法[J].計(jì)算機(jī)工程與應(yīng)用，2016，52（13）：152-156，182.

[4] 周利榮，胡天磊.Demytko素?cái)?shù)構(gòu)造算法優(yōu)化及應(yīng)用研究[J].電腦編程技巧與維護(hù)，2018（6）：54-59.

[5] Wilfrid Keller.Fermat numbers website data[EB/OL].[2020-12-20].http：//www.prothsearch.com/fermat.html.

[6] Ishmukhametov S T，Rubtsova R，Savelyev N.The error probability of the miller-Rabin primality test[J].Lobachevskii Journal of Mathematics，2018，39（7）：1010-1015.

[7] 魏成行.素性檢測算法研究及其在現(xiàn)代密碼學(xué)中的應(yīng)用[D].濟(jì)南：山東大學(xué)，2009.

【通聯(lián)編輯：謝媛媛】