求模逆元的幾種算法

摘要：基于模乘法逆元的定義、存在條件及其相關定理，首先，對各求模逆元的算法思想和計算過程進行了深入的剖析，并總結了它們各自的運算特點以及它們的局限性所在，最后，依據可計算的復雜性理論和實際所測試的數據，比較了各種算法的執行效率以及它們的使用范圍。

關健詞：模逆元；擴展歐幾里得算法；二進制擴展歐幾里得算法；牛頓迭代法；費馬小定理

中圖分類號：TP301文獻標識碼：A文章編號：1009-3044(2008)11-20308-03

1 引言

模算術就是用算術表達式模一些非零整數的計算。其數學概念就是是剩余類環。在數論和近世代數中，模算術是一個很重要的概念和方法，而求解模逆元又是模運算中最重要的一個問題。隨著計算機網絡技術和通信技術的發展，求解模逆元的問題在現代信息密碼安全中有著充分的應用和體現。尤其是在現代公鑰密碼體制中，加密和解密的過程中，都涉及到求逆元的問題。因此，研究求解模逆元的算法問題不僅在代數學上有重要的理論意義，而且在現代密碼學方面也有重要的現實意義。

2 模逆元概念及其相關定理

定義若整數a，b，p， 滿足a·b≡1(mod p)．則稱a為b模p的乘法逆元，即a=b-1mod p．其中，p是模數。

定理1設p是素數，a∈Z，則ap≡a mod p.特別地當gcd(a，p)=1時，

ap-1≡1(mod p)

成立。此定理就是有名的費馬小定理。作為代數學和數論學科方面的一個強有力的工具，它在多項式分解和素數檢測方面已有廣泛地應用。在此，進行適當轉換我們可以十分方便地利用它來計算mod p的逆元。即

ap≡a mod p=>a·ap-2≡1mod p.

所以， a-1=ap-2mod p.

因此，在求a mod p逆元時，若p為素數且gcd(a，p)=1時，則可直接利用該定理計算ap-2mod p的大小。其值就是所求的逆元a-1。一般情況下，由于p值較大，其計算煩瑣且易出錯。但由于計算思想相對簡單即先進行冪運算而后再模運算。因此，更一般的做法是：我們可以利用某些數學工具軟件(如maple)進行直接計算即可。

定理2 對于a∈{0，1，2，…，p-1}，存在a-1∈{1，2，3，…，p-1}使

a·a-1≡1mod p

充要條件是gcd(a，p)=1.

該定理的證明(此處從略)主要是用了最大公約數的有關性質，其中重要的一個就是，

gcd(a，p)=S·a+t·p(*1)

(*1)式表明，兩整數的最大公約數可用整系數線性組合的方式來表示。且當gcd(a，p)=1時，

有 s·a≡1mod p=>a-1=s.

該定理不但給出了一般模p逆元的存在條件，而且也提供了求逆元的一種方法。本文依據模逆元的定義和相關定理，下面依次給出了求模逆元的幾種算法。

3 求模逆元的算法

3.1 窮舉法求逆元

在a∈{1，2，…，p-1}的集合中，找某元素b，使得a·b≡1mod p成立的最直接方法：逐元素嘗試，直至找到符合條件的元素。此方法就是窮舉算法。其思想就是依據定義，在特定的集合中，逐一進行驗證。算法描述如下，

voidEnA(int a ，intp ，int inverse=0)// 窮舉法求逆元

{ for(int j←1;j

if (a*j % p = =1 ){ inverse←j; break;}

}//若inverse值是零，則表明a模p乘法逆元不存在。否則，inverse即為所求。

3.2 擴展歐幾里得方法求逆

歐幾里得算法是計算最大公約數的一個方便而又實用的方法。其原理是，重復利用帶余除法，直至余數為零時止。其過程可用下面的式子來描述。

gcd(p，a)=(a，r1)=…=gcd(rN-2，rN-1)=rN.

對上式，作適當處理：自后而前，逐次消去rN-1，rN-2，…，r2，r1

可得到一系列整數si，ti使得ri=sia+tip 且有rN=gcd(p，a)=sNa+tNp

上式就是前面提到的最大公約數的一個性質，即兩整數的最大公約數可用整系數線性組合的方式計算得到。這種計算相應的系數并通過線性組合的形式求最大公約數的方法，就是擴展歐幾里得算法。

關于sN，tN的計算可由下式遞歸得到（其正確性可用歸納法予以證明，此處從略），

其中，q表示每次輾轉相除時所得的商。顯然，擴展歐幾里得算法是在歐幾里得算法的基礎上，通過適當的變換處理得到。但比歐幾里得算法提供了更多的信息。如當gcd(a，p)=1時，a-1=sN.并且相應的系數都被一一記錄下來。相應算法描述如下，

int EEA(int a，int p) //擴展歐幾里得算法

{// 初始化

r0←a，r1←p，s0←1，s1←0，t0←0，t1←1，i←1;

while (ri≠0) //迭代計算

{qi←ri-1/ri；ri+1←ri-1%ri；si+1←si-1-qi·si；

ti+1←ti-1-qi·ti；i←i+1;

}

return .ri，si，ti；//其中si即為所求.

}

3.3 二進制擴展歐幾里得方法求逆

二進制擴展歐幾里得算法是在擴展歐幾里得算法的基礎上，基于以下性質作了改進。

(1)若a，b都是偶數，則有 gcd(a，b)=2gcd(a/2，b/2 )

(2)若a，是偶數，b是奇數，則有 gcd(a，b)=gcd(a／2，b )

(3)若a，b都是奇數，則有 gcd(a，b)=gcd((a-b)/2，b)

算法描述如下，

int BEEA (int a，int p)//二進制拓展算法

{ if (a==b) return (1，0);

if(both a and b are even) return EEA(a/2，b/2);

if(exactly one of the two，say a is even)

{(si，ti)←EEA(a/2，b);

if(si is even) return (si/s，ti);

else return (si+b)/2，ti-a/2);

}

if(both a and b are odd，say a>b)

{ (si，ti)←EEA((a-b)/2，b);

if (si is even) return (si/2，ti-si/2);

else return ((si+b)/2，ti-(si+a)/2);

}

}

該算法是擴展歐幾里得算法的一個變體。它主要是利用了公約數的有關性質并進行適當轉換及通過利用移位操作比用除法執行更快的特點對擴展歐幾里得算法進行了改進而得到的。

3.4 牛頓迭代法求逆

k=0，1，2，3，… （*2）

(*2)式就是求方程數值解的牛頓迭代公式。其計算過程是從合適的初值開始，依次代入公式就得到一個序列，該序列在迭代過程中會趨向一個穩定值，該值就是所求方程的解。其思想本質，就是用線性化方法去逼近方程的解。為此，我們可以把求a mod pl逆元看作是求方程a·x≡1mod pl的解。故，可通過找到1/x-a=0的一個近似根進行替代。將其代入(*2)得

（*3）

選擇合適的初值x0(滿足a·x≡1mod p)，利用(*3)，就可求形如a mod pl的逆元。算法描述如下：

int NIA(int a，int p，int l)//牛頓迭代法求逆元

{x0←EEA(a，p，temp);//x0:迭代初值

r←log2 l //迭代的次數

for(i←1;i≤r;i++){xi←(2xi-1-a·xi-1 2)modp2i;}

return xr ;

}

4 算法性能分析與比較

4.1 定性分析

窮舉法求逆的時間復雜度為O(p)；擴展算法的時間復雜度為O(log(max(a，p)))；二進制擴展算法的時間復雜度為O(log(max(a，p)))，但其中參與運算的常數因子與擴展歐幾里得算法中常數因子相比會更小；牛頓迭代法求逆的時間復雜度為O(l·log(p))。其中，p為模數，max(a，p)指取兩數之中的大數，l為p的指數次冪。根據復雜度的函數表達式，從理論上可得出結論：隨著p，a的增大，窮舉法時間消耗增長幅度最大，牛頓迭代法次之，擴展算法有緩慢增加，二進制擴展算法的時間消耗增加將最為緩慢。

4.2 定量比較

為了便于比較各求模逆元算法的時間開銷情況，通過實際測試，我們在下表中列出了相應的測試結果。具體操作環境是，WinXP操作系統，Visual C++6.0，Pentium(R)4 CPU 1.8GHz， 256MB內存。

各求模逆元算法的時間開銷對比表

注：（1）重復次數是指執行算法函數的次數；

（2）平均時間是指每次運行算法函數所消耗的時間。

5 結束語

理論分析和實際的數據測試表明：擴展歐幾里得算法是求逆元的基本算法，二進制擴展算法是最優的一種算法。在模數較小的情況下用窮舉法也是可行的。注意牛頓迭代法以求解方程的方法計算逆元的思想是相比于其它算法是個不錯的方法。但由于在計算過程中，含有較多的乘法和除法運算，致使其時間消耗較多。另外，根據實際情況利用費馬小定理(借助某些工具軟件)進行快速地求解和驗證也是一個很便利的方法。

參考文獻：

[1] J von zur Gathen，J Gerhard.Modern Computer Algebra[M].北京:世界圖書出版公司，2001.

[2] (美)Michael T Goodrich，RobertoTamassia.霍紅衛，譯.算法分析與設計[M].北京:人民郵電出版社，2006.

[3] (美)薩尼(Sahni S).汪詩林，等，譯.數據結構、算法與應用[M].北京:機械工業出版社，2000.

注：本文中所涉及到的圖表、注解、公式等內容請以PDF格式閱讀原文