一種SM9 算法R-ate 對的快速實現(xiàn)方法*

胡芯憶, 何德彪, 彭 聰, 羅 敏, 黃欣沂

1. 武漢大學 國家網絡安全學院 空天信息安全與可信計算教育部重點實驗室, 武漢 430072

2. 香港科技大學(廣州) 信息樞紐 人工智能學域, 廣州 511455

1 引言

基于身份的密碼體制(identity-based cryptography, IBC) 由Shamir[1]于1984 年提出, 其思想是將用戶的身份信息作為公鑰, 由可信的密鑰生成中心(key generation center, KGC) 為每個用戶生成其身份信息對應的私鑰, 這種方式無需使用公鑰證書來證明用戶身份與用戶公鑰的綁定關系, 簡化了公鑰基礎設施(public key infrastructure, PKI) 中用戶公鑰證書的管理. 2016 年3 月, 國家密碼管理局發(fā)布行業(yè)標準“GM/T 0044-2016《SM9 標識密碼算法》”, 包括總則、數(shù)字簽名算法、密鑰交換協(xié)議、密鑰封裝機制和公鑰加密算法、參數(shù)定義五個部分[2]. 成為國家標準以后, SM9 算法在各行業(yè)領域的應用更加廣泛. 雙線性對是SM9 密碼算法的核心運算部件, 它的性能決定了SM9 密碼算法的性能. 研究雙線性對的高效實現(xiàn)對推動SM9 算法的應用和落地具有重要的價值.

常見的雙線性映射主要有Weil 對[3]和Tate 對[4]. 雙線性對可以通過Miller 算法[5]有效計算.研究發(fā)現(xiàn)Tate 對在計算效率方面優(yōu)于Weil 對[6,7], 因此Tate 對一直是研究關注的重點. 后來, 學者們相繼提出了Tate 對的變體, 如Ate 對[8]、O-ate 對[9]和R-ate 對[10]等, 這些雙線性對有效降低了Miller 算法的迭代次數(shù), 減少了計算量, 計算效率更高. Barreto 等人[11]提出了Tate 對新的定義公式, 減少了Miller 算法中的無關運算, 避免了擴域中的求逆運算. Barreto 和Naehrig 提出了一族配對友好(pairing-friendly) 的素數(shù)階橢圓曲線, 嵌入度為12, 其定義為E:y2=x3+b,b/= 0, 這類曲線統(tǒng)稱為Barreto-Naehrig 曲線, 即BN 曲線[12], 可以有效降低Miller 算法的迭代次數(shù). Scott 等人[13]利用非超奇異橢圓曲線上的有效自同態(tài), 通過計算Miller 算法的中間值減少Tate 對的計算量. Galbraith等人[14]在大素數(shù)特征域上的一大類橢圓曲線上構造了一種可計算的有效自同態(tài), 可用于加速雙線性對中Frobenius 映射的計算. Costello 等人[15]減少了Ate 雙線性對Miller 算法的循環(huán)次數(shù), 并運用到高階扭曲線上進行運算.

冪指數(shù)運算(final exponentiation, FinalExp) 是雙線性對運算中另一個耗時的操作. Devegili 等人[16]優(yōu)化了BN 曲線上的FinalExp 運算并給出了Tate 對和Ate 對的快速實現(xiàn), 比之前的方法更快速、占用的內存更少. Scott 等人[17]通過分析多種配對友好曲線的結構, 提出了Tate 對及其變體的FinalExp 困難部分的通用實現(xiàn)方法, 比之前的方法更快速、更加通用. Granger 等人[18]提出了域特征滿足p ≡1 (mod 6) 的分圓子群中平方運算的快速實現(xiàn), 大大提高了FinalExp 困難部分的計算效率. Fuentes 等人[19]提出雙線性對指數(shù)運算后仍是雙線性對, 將FinalExp 困難部分的指數(shù)替換為另一個容易分解的指數(shù), 從而減少加法鏈的運算量. Duquesne 等人[20]針對BN 曲線上Tate 對及其變體的FinalExp 困難部分提出了4 種新的計算方法, 與以前的實現(xiàn)方法相比需要的內存資源減少了約37%.

Beuchat 等人[21]設計了BN 曲線上O-ate 對的快速軟件實現(xiàn)庫, 通過選擇BN 曲線的參數(shù)來減少Miller 算法和FinalExp 的操作數(shù)量, 在Intel Core i7 @2.8 GHz 處理器平臺, 計算254 比特BN 曲線下的O-ate 對需要0.832 毫秒. Naehrig 等人[22]提出了有限域元素新的表示方法, 利用AMD64 架構的雙精度浮點SIMD 指令在257 比特的BN 曲線上實現(xiàn)了O-ate 對, 在英特爾Core 2 Quad Q6600 處理器上計算一次O-ate 對需要4470 408 個時鐘周期. Aranha 等人[23]將延遲約減(lazy reduction) 應用到整個雙線性對計算過程, 提出當BN 曲線參數(shù)為負數(shù)時避免FinalExp 中求逆運算的方法, 在Phenom I 平臺計算BN254 曲線下的O-ate 對需要個1562 000 個時鐘周期. Azarderakhsh 等人[24]比較了使用仿射坐標和標準射影坐標實現(xiàn)O-ate 對的計算效率, 并在多個平臺實現(xiàn)了多個安全等級BN 曲線下的O-ate對, 在x86-64 AMD Opteron II @2.4 GHz 處理器平臺, 使用C 語言和標準射影坐標計算BN254 曲線下的O-ate 對需要3 290 000 個時鐘周期. 針對SM9 算法雙線性對, Zhen 等人[25]分析了SM9 雙線性對的計算過程和相應的優(yōu)化方法, 將Miller 循環(huán)和FinalExp 的效率分別提高了約5.2% 和0.91%. Cheng等人[26]基于SM9 算法參數(shù)對R-ate 雙線性對進行了優(yōu)化, 并通過Frobenius 映射來降低硬件成本, 基于C++ 驗證了優(yōu)化算法的有效性, 并基于FPGA 提出了硬件并行設計.

本文工作. 本文結合當前雙線性對的優(yōu)化實現(xiàn)方法, 提出了SM9 算法雙線性對的快速實現(xiàn)方法. 理論分析顯示, 本方法的計算復雜度低于現(xiàn)有的實現(xiàn)方法. 實驗結果顯示, 在Intel(R) Core(TM) i7-6500U@2.50 GHz 處理器上, 本方法計算一次SM9 算法雙線性對的時間為3.389 毫秒.

2 預備知識

2.1 符號定義

表1 給出了文中使用符號的定義.

表1 符號定義Table 1 Notations and their definitions

2.2 BN 曲線

根據(jù)文獻[11], BN 曲線E在有限域Fp上點群E(Fp) 的階為r= #E(Fp),r是一個素數(shù). 基域特征p, 群的階r和Frobenius 跡tr都由參數(shù)t定義:

其中t可以是滿足p(t) 和r(t) 都為素數(shù)的任意整數(shù). 根據(jù)《SM9 標識密碼算法參數(shù)定義》, SM9 算法使用BN 曲線, 參數(shù)t=0x600000000058F98A,p和r為256 比特長的素數(shù). 由于BN 曲線的嵌入度為12,SM9 算法雙線性對需要在有限擴域Fp12中進行計算.

2.3 R-ate 對

算法1 R-ate 對),s = 6t+2 Output: Ra(Q,P)1 設s = ∑L-1 i=0 si2i, sL-1 = 1, si ∈{0,1};Input: P ∈E(Fp), Q ∈E′(Fp2 2 初始化T ←Q,f ←1;3 for i from L-2 to 0 do 4 計算f ←f2 ·gT,T (P); T ←[2]T;5 if si = 1 then f ←f ·gT,Q(P); T ←T +Q;7 end 8 end 9 計算Q1 ←πp(Q), Q2 ←πp2 (Q);10 計算f ←f ·gT,Q1 (P); T ←T +Q1;11 計算f ←f ·gT,-Q2 (P); T ←T -Q2;12 計算f ←fp12-1 6 r ;13 輸出f.

3 擴域Fp12 的優(yōu)化平方算法

3.1 塔式擴張構建方法

塔式擴張是R-ate 對運算的基礎, 根據(jù)《SM9 標識密碼算法參數(shù)定義》[2], Fp12采用1-2-4-12 的塔式擴張方式構建, 即:

定理1 1-2-4-12 擴張方式和1-2-6-12 擴張方式可以相互轉換.

證明: 設f是擴域Fp12的元素. 根據(jù)1-2-4-12 的擴張方式, 假設f=a+bw+cw2=(a0+a1v)+(b0+b1v)w+(c0+c1v)w2=a0+b0w+c0w2+a1w3+b1w4+c1w5. 根據(jù)1-2-6-12 的擴張方式, 假設f=g+ht= (g0+g1s+g2s2)+(h0+h1s+h2s2)t=g0+h0t+g1t2+h1t3+g2t4+h2t5. 根據(jù)兩種擴張方式的二次擴域構建是相同的, 有w6=t6=u, 且a0+b0w+c0w2+a1w3+b1w4+c1w5=g0+h0t+g1t2+h1t3+g2t4+h2t5, 因此有各項系數(shù)對應相等, 得到系數(shù)轉換如圖1 所示. 同理, 也能得到1-2-6-12 擴張到1-2-4-12 擴張的轉換.

圖1 1-2-4-12 擴張到1-2-6-12 擴張的轉換Figure 1 Conversion from 1-2-4-12 to 1-2-6-12 extension

可以發(fā)現(xiàn), 兩種擴張方式的轉化是沒有開銷的, 因此我們可以利用擴域Fp6中的快速平方和求逆來加速最后冪運算中分圓域Gφ6(Fp2) 的運算, 文獻[18] 給出了具體的算法.

3.2 優(yōu)化平方算法

本文采用文獻 [18] 的方法來實現(xiàn)擴域 Fp2的平方運算, 假設a=a0+a1u ∈Fp2, 有a2= (a0+βa1)(a0+a1)-a0a1(1 +β) + 2a0a1u. 類似地, 該方法也可用于擴域Fp4的平方運算. 本文將文獻[27] 的SQR3 方法用于擴域Fp12元素的平方運算, 如算法2 所示, 該算法需要4 次Fp4平方、1 次Fp4乘法和一些加法運算. 表2 總結了基域及各有限擴域中基礎運算的計算開銷, 在實現(xiàn)Fp12平方運算時, 本文的方法需要11 個Fp2域上的乘法, 比文獻[21,23,24] 的方法少一個.

算法2 擴域Fp12 的平方運算Input: f = f0 +f1w+f2w2 ∈Fp12 Output: h = f2 = h0 +h1w+h2w2, 其中h0 = f20 +2f1f2v, h1 = f22 v+2f0f1, h2 = f21 +2f0f2 1 h0 ←f20;2 h1 ←f22;3 s0 ←f2 +f0;4 s1 ←(s0 -f1)2;5 s0 ←(s0 +f1)2;6 s2 ←f1f2;7 s2 ←2s2;8 s3 ← s0+s1 2 ;9 h2 ←s3 -h1 -h0;10 h1 ←h1v+s0 -s2 -s3;11 h0 ←h0 +s2v;12 輸出h = h0 +h1w+h2w2.

表2 基域和有限擴域的計算開銷Table 2 Computational cost of based field and finite extension fields

4 稀疏乘法

當x1/=x2時, 設gT,Q是穿過T和Q的直線, 則gT,Q=y-mx+mx1-y1. 因此gT,Q在點P(xP,yP)∈E(Fp):y2=x3+b(b/=0) 的值為

當x1=x2時, 類似地, 設gT,T是過T點的切線, 其在點P(xP,yP) 的值為

對于六次乘扭曲線, 線函數(shù)在P點的值為g= ?g0+ ?g2w2+ ?g3w3稀疏元素的形式. 對于f=f0+f1w+f2w2,fi ∈Fp4,g= ?g0+?g2w2+?g3w3=g0+?g2w2,?gi ∈Fp2, 設h=f ·g, 則

算法3 給出了乘扭曲線下用稀疏乘法計算f ·g的算法描述.

算法3 乘扭曲線稀疏乘法Input: f = f0 +f1w+f2w2,g = ?g0 + ?g2w2 + ?g3w3 = g0 + ?g2w2 Output: h = f ·g = h0 +h1w+h2w2 1 t0 ←f0 ·g0;2 t1 ←Fp4SparseMul(f2,?g2);3 u0 ←Fp4SparseMul(f1 +f2,?g2);4 u1 ←(f0 +f2)·[(?g0 + ?g2)+ ?g3v];5 u2 ←(f0 +f1)·g0;6 h0 ←t0 +(u0 -t1)v;7 h1 ←u2 -t0 +t1v;8 h2 ←u1 -t0 -t1;9 輸出h = h0 +h1w+h2w2.

類似地, 對于六次除扭曲線, 線函數(shù)在P點的值為g= ?g0+?g1w+?g3w3稀疏元素的形式[24]. 對于f=f0+f1w+f2w2,fi ∈Fp4,g= ?g0+?g1w+?g3w3=g0+?g1w,?gi ∈Fp2, 設h=f ·g, 則

算法4 給出了除扭曲線下用稀疏乘法計算f ·g的算法描述.

算法4 除扭曲線稀疏乘法Input: f = f0 +f1w+f2w2,g = ?g0 + ?g1w+ ?g3w3 = g0 + ?g1w Output: h = f ·g = h0 +h1w+h2w2 1 t0 ←f1 ·g0;2 t1 ←Fp4SparseMul(f2,?g1);3 u0 ←Fp4SparseMul(f0 +f2,?g1);4 u1 ←(f1 +f2)·[(?g0 + ?g1)+ ?g3v];5 u2 ←(f0 +f1)·g0;6 h0 ←u2 -t0 +t1v;7 h1 ←u0 +t0 -t1;8 h2 ←u1 -t0 -t1;9 輸出h = h0 +h1w+h2w2.

算法3 和算法4 都調用了函數(shù)Fp4SparseMul, 將一個稀疏的Fp4元素(即Fp2元素) 與一個Fp4元素相乘, 具體如算法5 所示. 算法3 和4 需執(zhí)行3 次Fp4乘法、2 次Fp4SparseMul、9 次Fp4加法和1 次Fp2加法,根據(jù)表2,一次稀疏乘法的計算開銷為3·(3 ?m+5?a+mξ)+2·2 ?m+9·2?a+?a=13 ?m+34?a+3mξ.

算法5 Fp4SparseMul Input: a = a0 +a1v ∈Fp4,b0 ∈Fp2 Output: c = a·b0 = c0 +c1v 1 c0 ←a0 ·b0;2 c1 ←a1 ·b0;3 輸出c = c0 +c1v.

5 Miller 算法優(yōu)化實現(xiàn)

5.1 基于NAF 的點乘運算優(yōu)化實現(xiàn)

因此對算法1 的第3 步改進如下, 將s表示為2-NAF 形式, 從次高位至低位掃描每一位si. 每次計算點倍T= [2]T, 求線函數(shù)值, 并計算f2·gT,T(P); 若si= 1, 則計算點加T+Q, 求線函數(shù)值, 并計算f ·gT,Q(P); 若si=ˉ1, 則計算點加T+(-Q), 求線函數(shù)值, 并計算f ·gT,-Q(P). 為此, 我們需要預先計算-Q, 該操作只需要Fp2上的取負. 改進的雙線性計算如算法6 所示. 原始的s比特長度為66, 漢明重量為16, 2-NAF 形式的s比特長度為66, 漢明重量為11, 因此將s表示為2-NAF 形式可以減少5 次點加運算、線函數(shù)求值和稀疏乘法.

5.2 Frobenius 映射優(yōu)化實現(xiàn)

算法6 改進的R-ate 對Input: P ∈E(Fp), Q ∈E′(Fp2),s = 6t+2 Output: Ra(Q,P)1 設s = ∑L-1 i=0 si2i, sL-1 = 1, si ∈{ˉ1,0,1};2 初始化T ←Q,f ←1;3 for i from L-2 to 0 do f ←f ·gT,Q(P); T ←T +Q;7 end 8 else if si = ˉ1 then 4 計算f ←f2 ·gT,T (P); T ←[2]T;5 if si = 1 then 6 f ←f ·gT,-Q(P); T ←T -Q;10 end 11 end 12 計算Q1 ←πp(Q), Q2 ←πp2 (Q);13 計算f ←f ·gT,Q1 (P); T ←T +Q1;14 計算f ←f ·gT,-Q2 (P); T ←T -Q2;15 計算f ←fp12-1 9 r ;16 輸出f.

5.3 Miller 算法計算復雜度分析

根據(jù)文獻[24], 我們使用標準射影坐標計算點運算和線函數(shù), 使用該坐標系, 一次點倍及線函數(shù)求值的計算開銷為2 ?m+7?s+22?a+4m, 一次點加及線函數(shù)求值的計算開銷為11 ?m+2?s+8?a+4m. 通過以上的優(yōu)化, Miller 算法中一共需要執(zhí)行65 次Fp12平方、點倍及線函數(shù)求值、稀疏乘法, 12 次點加及線函數(shù)求值、稀疏乘法,2 次Frobenius 映射. 根據(jù)表2,Miller 算法的計算開銷為65·[(11 ?m+43?a+11mξ)+(2 ?m+7?s+22?a+4m)+(13 ?m+34?a+3mξ)]+12·[(11 ?m+2?s+8?a+4m)+(13 ?m+34?a+3mξ)]+2·(2a+4m)=1978 ?m+479?s+6939?a+312m+946mξ+4a.

6 冪指數(shù)優(yōu)化實現(xiàn)

根據(jù)文獻[21,29],FinalExp 中3 次指數(shù)為t的模冪運算約占FinalExp 運算時間的79%,而FinalExp約占整個雙線性對計算時間的42%. 因此, 優(yōu)化冪指數(shù)運算可以有效提高R-ate 對的性能.

6.1 基于NAF 的冪指數(shù)優(yōu)化實現(xiàn)

2-NAF 形式的t比特長度為63, 漢明重量減少為11, 從而在一次計算mt時, 可以減少3 次擴域Fp12乘法運算.

我們采用文獻[29] 提出的壓縮平方算法來提高平方運算的效率. 對于g=(g0+g1v)+(g2+g3v)w+(g4+g5v)w2∈Gφ6(Fp2), 壓縮函數(shù)C和解壓縮函數(shù)D定義如下[29]:

其中,

若h=g2=(h0+h1v)+(h2+h3v)w+(h4+h5v)w2, 那么有D(g2)=[h2,h3,h4,h5], 根據(jù)文獻[29] 的SQR2345, 其中

根據(jù)Karatsuba 乘法[30], 6g4g5=3·[(g4+g5)2-g24-g25], 6g2g3=3·[(g2+g3)2-g22-g23]. 因此解壓縮函數(shù)的計算開銷是3?s+2 ?m+?i+9?a+2mξ, 壓縮乘法的計算開銷是6?s+28?a+3mξ.

改進的模冪運算如算法7 所示, 利用Montgomery 同時求逆方法[24], 指數(shù)為t的冪指數(shù)運算的計算開銷為62·(6?s+28?a+3mξ)+10·(18 ?m+60?a+8mξ)+11·(3 ?m+3?s+9?a+2mξ)+3·10 ?m+?i= 243 ?m+ 405?s+?i+ 2435?a+ 288mξ. 若不采用壓縮平方, 指數(shù)為t的冪指數(shù)運算的計算開銷為62·(6 ?m+37?a+7mξ)+10·(18 ?m+60?a+8mξ)+3?a=552 ?m+2897?a+514mξ, 因此采用壓縮平方算法能夠減少約6.8% Fp上的乘法和約16% Fp上的加法運算.

算法7 基于2-NAF 的模冪運算Input: m,m-1,t = ∑l-1 i=0 ti2i,ti ∈{ˉ1,0,1},tl-1 = 1 Output: c = mt 1 初始化c ←t;2 for i from L-2 to 0 do 3 計算c ←c2;4 if si = 1 then c ←c·m;6 end 7 else if si = ˉ1 then 5 c ←c·m-1;9 end 10 end 11 輸出c.8

6.2 FinalExp 計算復雜度分析

對于簡單部分,fp6-1= ˉf · f-1[21], 因此需要一次共軛、一次求逆和一次乘法, 計算開銷為3?a+(25 ?m+9?s+?i+61?a+12mξ)+(18 ?m+60?a+8mξ)=43 ?m+9?s+?i+124?a+20mξ; 計算fp2+1, 需要一次p2次的Frobenius 映射和一次乘法, 計算開銷為5·2m+(18 ?m+60?a+8mξ)=18 ?m+60?a+10m+8mξ.

通過以上優(yōu)化, 對于困難部分, 需要計算ft,ft2,ft3,fp,fp2,fp3,(ft)p,(ft2)p,(ft3)p,(ft2)p2, 這10個值的計算開銷是3 次指數(shù)為t的模冪運算、4 次p次Frobenius 映射和3 次p2次Frobenius 映射,即3·(243 ?m+405?s+?i+2435?a+288mξ)+4·(10m+6a)+3·10m= 729 ?m+1215?s+3?i+7305?a+70m+864mξ+24a. 然后我們使用文獻[17] 的方法, 即加法鏈, 來計算FinalExp 的困難部分, 計算開銷為258 ?m+943?a+132mξ. 困難部分的計算開銷為987 ?m+1215?s+3?i+8248?a+70m+996mξ+24a.

因此, FinalExp 的計算開銷為1048 ?m+1224?s+4?i+8432?a+80m+1024mξ+24a.

7 實驗與結果

本節(jié)針對所提出的R-ate 對快速實現(xiàn)方法進行性能評估. 根據(jù)前文的分析, 表3 中總結了R-ate 對所需的計算開銷. 進一步地, 為了直觀地說明計算R-ate 對需要的運算數(shù)量, 我們將表3 中擴域Fp2的運算用基域Fp的乘法和加法運算表示. 依據(jù)《SM9 標識密碼算法參數(shù)定義》, 本文使用R-ate 對作為SM9算法雙線性對的實例. 由于R-ate 對的實現(xiàn)文獻較少, 考慮到O-ate 對與R-ate 對算法類似, 我們將本文的優(yōu)化實現(xiàn)方法用于O-ate 對, 并與O-ate 對的經典實現(xiàn)方法進行對比, 結果如表4 所示.

表3 R-ate 對所需的運算數(shù)量Table 3 Operation counts for R-ate pairing computation

表4 雙線性對實現(xiàn)所需的運算數(shù)量比較Table 4 Comparison of operation counts for pairing implementation

實現(xiàn)O-ate 對時, Miller 算法所需的Fp乘法數(shù)量比文獻[21,23,24] 更少, 比文獻[21] 和[24] 分別減少了約6.9% 和4.6%; FinalExp 所需的乘法數(shù)量比文獻[21] 減少了約18.6%, 但比文獻[23,24] 多, 這是由于文獻[23,24] 使用了文獻[19] 的技巧來優(yōu)化FinalExp, 使用另一個易于分解的指數(shù)來簡化加法鏈. 文獻[19] 保持了雙線性對的雙線性性質, 但改變了其運算結果. 本文嚴格遵循《SM9 標識密碼算法》[2]中雙線性對的參數(shù)選取, 為保證計算結果與標準示例的一致性, 未采用該方法進行優(yōu)化.

本文測試程序的構建基于密碼開源庫RELIC 0.5.0[31], 依據(jù)《SM9 標識密碼算法參數(shù)定義》使用C語言進行實現(xiàn). 測試環(huán)境為一臺配置有64 位Windows 10 操作系統(tǒng), Intel(R) Core(TM) i7-6500U CPU@2.50 GHz, 8.00 GB RAM 的惠普筆記本電腦. 取500 次實驗運行時間的平均值作為測試結果, 如表5 所示, 結果顯示本實現(xiàn)計算一次R-ate 對需要的時間為3.389 毫秒.

表5 R-ate 對實現(xiàn)性能測試Table 5 Test on efficiency of R-ate pairing

8 結論

雙線性對的計算效率對實現(xiàn)基于雙線性對的密碼協(xié)議具有非常重要的影響, 本文通過詳細分析SM9算法R-ate 對的Miller 算法和冪指數(shù)運算中的耗時操作, 提出了R-ate 對的一系列優(yōu)化算法, 包括優(yōu)化的十二次擴域平方算法、稀疏乘法、基于NAF 的點乘運算優(yōu)化實現(xiàn)、基于NAF 和加法鏈的冪指數(shù)優(yōu)化實現(xiàn). 通過與現(xiàn)有實現(xiàn)方法的對比分析, 本方法在計算復雜度方面具有較高的優(yōu)勢. 實驗結果顯示, 本實現(xiàn)計算一次R-ate 對需要的時間為3.389 毫秒.