整數(shù)上離散高斯取樣的常數(shù)時間實現(xiàn)方法

盧嘉嘉,杜育松

(中山大學(xué) 數(shù)據(jù)科學(xué)與計算機(jī)學(xué)院,廣州 510006)

0 概述

格密碼被認(rèn)為是未來能夠抵抗量子計算機(jī)攻擊的一類新型密碼[1-2]。研究表明,使用量子計算機(jī)有可能破解基于大整數(shù)分解或離散對數(shù)問題的傳統(tǒng)公鑰密碼,如RSA和橢圓曲線公鑰密碼。常見的數(shù)字簽名和公鑰加密等密碼學(xué)原語,都已經(jīng)有了基于格密碼的實例[3-4]。

離散高斯取樣是指在n維格中取出格點,滿足特定參數(shù)的n維離散高斯分布[5]。離散高斯取樣在格密碼中具有基礎(chǔ)性作用,是許多格密碼體制實現(xiàn)的基本操作[6-7]和決定這些密碼體制安全性的重要因素[8]。整數(shù)上的離散高斯取樣是一般離散高斯取樣的子問題,也是實現(xiàn)一般離散高斯取樣的核心子過程[5,9],主要指產(chǎn)生符合含特定參數(shù)的一維離散高斯分布的整數(shù)[10]。由于整數(shù)上的離散高斯取樣相比于一般離散高斯取樣,具有計算簡單且取樣效率高的特點,因此在一些密碼體制中直接作為一般n維格上的離散高斯取樣被使用[11]。然而,傳統(tǒng)的連續(xù)高斯分布取樣算法不能直接應(yīng)用于這種離散的情形[12],現(xiàn)有的整數(shù)上的離散高斯取樣算法預(yù)計算存儲量大,無法抵御計時攻擊,或者取樣性能不高,無法很好地滿足格密碼實現(xiàn)的要求[13],一定程度上制約了格密碼方案的推廣與應(yīng)用。因此,研究整數(shù)上的離散高斯取樣算法對于格密碼的實現(xiàn)具有重要的理論和實際意義。

針對密碼方案的側(cè)信道攻擊是密碼學(xué)領(lǐng)域非常重要的研究課題。但格密碼體制有可能因為針對離散高斯取樣過程的計時攻擊導(dǎo)致秘密信息泄露[14-15],而先前設(shè)計的整數(shù)上的離散高斯取樣算法大多無法抵御計時攻擊。因此,在保證離散高斯取樣功能和效率的同時,如何防止離散高斯取樣由于計時攻擊造成的秘密信息泄露成為一個倍受關(guān)注的問題[16]。

基于Knuth-Yao算法的離散高斯取樣算法的優(yōu)點是取樣速度快,不需要在線的浮點運算,并且能夠?qū)崿F(xiàn)常數(shù)時間取樣。例如,文獻(xiàn)[17-18]將基于Knuth-Yao算法的取樣算法轉(zhuǎn)化為常數(shù)時間的對隨機(jī)比特串的位操作,既保證了取樣速度,又防止了由于計時攻擊造成的秘密信息泄露。然而,在實現(xiàn)該算法的過程中,將取樣過程轉(zhuǎn)化為對隨機(jī)比特串的位操作包含了復(fù)雜的預(yù)計算過程(需要編寫專門的復(fù)雜預(yù)計算程序去實現(xiàn),預(yù)計算中涉及的部分參數(shù)還需要通過額外的實驗數(shù)據(jù)比對才能確定),并且實現(xiàn)位操作的代碼規(guī)模巨大,有可能導(dǎo)致硬件實現(xiàn)中所需的電路規(guī)模較大,不適合在存儲資源受限的電子設(shè)備上實現(xiàn)。

針對整數(shù)上離散高斯取樣的問題,本文提出一種常數(shù)運行時間的簡單實現(xiàn)方法。該方法在支持單指令多數(shù)據(jù)(SIMD)的設(shè)備上采取向量化操作手段,從而提高取樣速度。

1 預(yù)備知識

1.1 離散高斯分布及取樣

1.2 Knuth-Yao算法

實際上在采樣操作之前不需要將DDG樹按樹的形式存儲,它可以從概率矩陣中被快速地構(gòu)造并使用。考慮樣本點概率的二進(jìn)制擴(kuò)展,它們可以寫為二進(jìn)制矩陣的形式,寫成的矩陣被稱為概率矩陣并由M表示。對于i,j≥0,設(shè)pj=pj0+pj12-1+pj22-2+…+pji2-i+…是第j個樣本點概率的二進(jìn)制展開,其中pji∈{0,1}。在概率矩陣M中,第j行Mj即為(pj0,pj1,…,pji,…),對應(yīng)第j個樣本點概率的二進(jìn)制展開。概率矩陣與DDG樹的對應(yīng)關(guān)系如圖1所示。

圖1 概率矩陣與對應(yīng)的DDG樹Fig.1 Probability matrix and corresponding DDG tree

基于對DDG樹內(nèi)部節(jié)點相對距離的分析,文獻(xiàn)[12]將基于Knuth-Yao采樣過程描述為如算法1所示的形式。

算法1基于離散高斯分布的Knuth-Yao取樣算法[4]

輸入λ+1行的概率矩陣M

輸出采樣結(jié)果x

1.j←0,d←0

2.采樣均勻隨機(jī)比特r

3.d←2d+r

4.for i←λ down to 0 do

5.d←d-M[i][j]

6.return i if d=-1,否則continue

7.end for

8.j←j+1并返回步驟2

在算法1中,d表示DDG樹訪問節(jié)點和同層最右節(jié)點的距離,j表示矩陣列索引,r取1表示往靠近內(nèi)部節(jié)點的方向游走。

對于一個具有無限支撐的離散分布,樣本點的數(shù)量及其概率的二元展開式可能是無限的,這將導(dǎo)致生成一個無限規(guī)模的概率矩陣。在現(xiàn)實中必須對概率矩陣進(jìn)行截斷處理,同時保證足夠的精度。如算法1中概率矩陣M的行數(shù)被設(shè)置為λ+1,使得算法輸出絕對值大于λ/2的概率是可忽略的。同時,矩陣M的列數(shù)也是有限的,因為每個樣本點的概率值只能以有限的精度存儲。一般來說,對于給定的σ,取λ=12σ(表示算法輸出的整數(shù)在-6σ～6σ之間),而矩陣的列數(shù)可以被設(shè)置為64 bit,從而保證大多數(shù)應(yīng)用的精度需求[10,19]。

1.3 窮盡比特掃描的Knuth-Yao算法

為抵御計時攻擊,可以將標(biāo)準(zhǔn)的Knuth-Yao取樣算法改造成為常數(shù)運行時間的“窮盡比特掃描”的Knuth-Yao取樣算法。具體地,要求執(zhí)行算法在得到d=-1時,并不立即返回輸出值,而是設(shè)法記錄輸出值,并繼續(xù)完成對整個概率矩陣的遍歷。

2 常數(shù)運行時間的實現(xiàn)方法

為減少運行時間,本文考慮對算法1進(jìn)行改進(jìn)。主要思想可以分為3點:快速初始化變量d;快速確定d=-1時所在的概率矩陣列,即d=-1時變量j的值;按比特掃描概率矩陣的第j列,返回輸出值。為實現(xiàn)上述方法,在預(yù)計算階段,需要確定并記錄概率矩陣各列的漢明重量(記為wt),從而借助概率矩陣各列的wt實現(xiàn)d的快速初始化和更新。例如,參數(shù)σ取6,則概率矩陣大小為73行64列。一列的wt較大概率不會超過64,因此各列的wt可以用6位二進(jìn)制數(shù)來表示。

2.1 快速初始化變量

在算法1中,變量d被初始化為0,以d←2d+r的方式進(jìn)行更新,然后按比特掃描概率矩陣,并在每掃描完一列時再以d←2d+r的方式進(jìn)行更新。

注意到一個給定概率矩陣的前若干列有可能均為零列,即wt=0。例如,當(dāng)參數(shù)σ=6且c=0時,概率矩陣的前4列均為零列。一般地,如果概率矩陣的前N列均為零列,那么對于j=0,1,…,N-1,算法1的第4步～第7步對變量d沒有影響。或者說,在掃描概率矩陣的前N列時,變量d僅以d←2d+r的方式進(jìn)行更新,且一定有d≥0。設(shè)執(zhí)行算法時產(chǎn)生隨機(jī)比特串為b0b1…bN-1…?？梢则炞C,在即將掃描概率矩陣的第N列時,變量d的值滿足:

d=b0×2N-1+b1×2N-2+…+bN-1×20

即變量d的值為:將比特串b0b1b2…bN-1看成是二進(jìn)制數(shù)時所表示的十進(jìn)制整數(shù)。于是,在概率矩陣的前N列均為零列時,可以直接將變量d初始化為(b0b1…bN-1)2,并令j=N,避免遍歷零列造成不必要的時間開銷。

2.2 d=-1時所在列的快速確定

算法1通過對概率矩陣的列逐比特遍歷以更新d的值。如果記錄了的概率矩陣每一列的wt值,從首列非零列開始,逐列檢查d-wt是否小于等于-1,當(dāng)首次發(fā)現(xiàn)d-wt≤-1時,就可以確定對應(yīng)的列即為輸出結(jié)果所在的列。相反,如果d-wt>-1,則令d←d-wt,再以d←2d+r的方式更新變量d的值,然后對概率矩陣下一列繼續(xù)檢查d-wt,直到在某一列時,發(fā)現(xiàn)d-wt≤-1。

稱輸出結(jié)果所在的列為目標(biāo)列,在確定目標(biāo)列的過程中,若找到滿足條件的列隨即停止掃描,攻擊者將有可能檢測到程序停止運行,通過計時攻擊的手段確定目標(biāo)列所在范圍,而目標(biāo)列的信息泄露進(jìn)而將導(dǎo)致輸出結(jié)果的信息泄露。為了防止這一問題,需在確定目標(biāo)列時,不立即停止檢查d-wt,而是設(shè)法記錄目標(biāo)列,繼續(xù)完成對剩余所有列的檢查,從而保證常數(shù)的運行時間。對上述確定目標(biāo)列的方法在下文以C/C++偽代碼的形式進(jìn)行描述。

1.WhichCol←FirstNonZeroCol,s←1

2.for j←FirstNonZeroCol to 63 do

3.d←d+s (d+(rnd64bits&1))

4.rnd64bits←rnd64bits?1

5.s←s&unsigned(d>=ColWts[j])

6.d←d-s ColWts[j]

7.WhichCol←WhichCol+s

8.end for

在偽代碼中,rnd64bits為預(yù)先生成的64位隨機(jī)比特串。FirstNonZeroCol表示概率矩陣的首列非零列。數(shù)組ColWts存儲概率矩陣各列的wt。引入變量s作為指標(biāo),并用WhichCol表示目標(biāo)列的位置。在確定目標(biāo)列的過程中,j從概率矩陣的首列非零列FirstNonZeroCol開始,直到j(luò)=63,檢查完概率矩陣后續(xù)所有列。

當(dāng)首次出現(xiàn)d

2.3 按比特掃描概率矩陣的第WhichCol列

在確定目標(biāo)列WhichCol和對應(yīng)的d值后,按比特掃描概率矩陣的第WhichCol列確定d=-1時所在的行WhichRow,即可確定返回的輸出值。為了確保常數(shù)運行時間,同樣引入變量s作為指標(biāo)。初始化s為1,在確定WhichRow后,s由1置為0,程序保持對第WhichCol列所有剩余比特的掃描。變量WhichCol雖然繼續(xù)參與運算,但它的值不再變化。

需要指出的是,為了抵御可能的cache攻擊[14],概率矩陣應(yīng)按行來存儲,而不應(yīng)按列來存儲。例如,參數(shù)σ取6時,概率矩陣有73行64列,這時應(yīng)使用73個uint64_t型整數(shù),而非64個bitset<73>("…")型對象來存儲整個概率矩陣。如果按列來存儲,因為程序會只讀取其中一個bitset<73>型對象,cache攻擊有可能確定當(dāng)前正在掃描的概率矩陣的第WhichCol列,從而造成輸出值的信息泄露。當(dāng)概率矩陣按行來存儲時,掃描概率矩陣的第WhichCol列需要概率矩陣所有行的數(shù)據(jù)參與計算,從而避免了cache攻擊的可能。

3 本文方法的實現(xiàn)與比較

在處理器i7-8550U內(nèi)存為16 GB的筆記本電腦上,使用g++編譯器進(jìn)行編譯,實現(xiàn)了本文設(shè)計提出的取樣過程,涉及的均勻隨機(jī)數(shù)通過AES-256-CTR方式產(chǎn)生。圖2比較了標(biāo)準(zhǔn)實現(xiàn)方法和本文實現(xiàn)方法的取樣速度。

圖2 參數(shù)σ控制下的算法采樣速度對比Fig.2 Comparison of algorithm sampling speed under parameter σ control

從圖2可以看出,本文的常數(shù)運行時間實現(xiàn)方法較算法1的標(biāo)準(zhǔn)實現(xiàn)方法速度略快,可以滿足密碼方案的性能需求。對于“窮盡比特掃描”版本的取樣算法,當(dāng)σ=6時,其取樣速度僅為0.155 5×106samples/s,并且隨著σ的增大,取樣速度將低至無法使用。

此外需要說明的是,圖2涉及的算法1的兩種實現(xiàn)方法在初始化變量d時,均使用了2.1節(jié)中的快速初始化方法,以避免遍歷零列造成的不必要的時間開銷。

本文算法中的概率矩陣行數(shù)被設(shè)置為12σ+1,實際上,對于中心化離散高斯分布,由于概率分布的對稱性,概率矩陣只需6σ+1行,這使得可以進(jìn)一步提高取樣算法的實現(xiàn)效率。本文為了保證普適性,考慮實現(xiàn)的是非中心化離散高斯分布的取樣。

在所有具有常數(shù)運行時間的、中心化離散高斯分布的取樣算法中,到目前為止,文獻(xiàn)[17]設(shè)計實現(xiàn)的取樣算法是速度最快的算法,達(dá)到15×106samples/s。本文算法的取樣速度雖然未能達(dá)到文獻(xiàn)[17]算法的速度,但是優(yōu)勢在于預(yù)計算過程簡單,存儲量小,對于資源受限的設(shè)備更加友好。

具體來說,本文的預(yù)計算過程僅是計算給定離散高斯分布的概率矩陣,并確定概率矩陣每個列向量的漢明重量。預(yù)計算存儲則只包括概率矩陣和每個列向量漢明重量的值。例如,對于參數(shù)為σ=6的中心化離散高斯分布,概率矩陣規(guī)模可以設(shè)置為73行64列,可以保證多數(shù)應(yīng)用所需的安全精度。73個行向量用73個uint64_t型整數(shù)存儲,而每個列向量漢明重量的值則可以用64個uint8_t型整數(shù)存儲,再用一個uint8_t型整數(shù)記錄一個非零列向量的位置,因此總的存儲量為73×8+64+1=649 Byte。如果需要更大的安全精度,則可以增加概率矩陣規(guī)模,如將概率矩陣的行和列數(shù)量都加倍,總的存儲量不超過4×649≈2.5 KB。

文獻(xiàn)[17]算法實現(xiàn)之前需要先將基于Knuth-Yao的取樣算法轉(zhuǎn)化為對隨機(jī)比特串的位操作(布爾表達(dá)式)。這一過程需要編寫復(fù)雜的預(yù)計算程序去實現(xiàn),涉及Karnaugh圖優(yōu)化技術(shù)或者Espresso啟發(fā)式邏輯最小化工具,并且預(yù)計算中涉及的部分參數(shù)需要通過額外的實驗數(shù)據(jù)比對才能確定。當(dāng)離散高斯分布的參數(shù)σ發(fā)生變化時,上述的整個轉(zhuǎn)化過程還需要重新進(jìn)行。因此,其整體預(yù)計算開銷較大。該算法在實現(xiàn)時,需要在源代碼中寫入大約18×103個位與運算 AND、5×103個位或運算 OR以及11×103個位非運輸 NOT[17]。如果每個運算占用一個字節(jié),那么源代碼規(guī)模將至少達(dá)到18+5+11=34 KB,遠(yuǎn)大于本文的2.5 KB。

4 算法的向量化實現(xiàn)

目前的CPU多數(shù)支持單指令多數(shù)據(jù)(SIMD),可以并行處理包含多個數(shù)據(jù)的向量。本文提出的取樣算法也可以通過使用SIMD將其進(jìn)行向量化,從而實現(xiàn)取樣速度的進(jìn)一步提高。特別地,可以使用Agner Fog的VCL (C++向量類庫)對取樣算法進(jìn)行向量化[20]。VCL向量化類庫在SIMD的支持下為整數(shù)和浮點數(shù)提供了優(yōu)化的向量操作,使程序員更容易編寫向量化代碼,而不必為調(diào)用SIMD指令而使用匯編語言或內(nèi)部函數(shù)。

VCL類庫支持Vec64c數(shù)據(jù)類型。如果AVX512指令集在目標(biāo)CPU上是可以使用的,那么Vec64c是具有64個int8_t型整數(shù)的向量類型。然而,目前大多數(shù)CPU支持AVX和AVX2指令集,但不支持AVX512。因此,本文僅使用了Vec32c類,它在AVX2的支持下并行處理具有32個int8_t型整數(shù)的向量。結(jié)合VCL類庫,可以將實現(xiàn)過程分為PopulatePoolColumns()和Sampler()兩個子函數(shù)。PopulatePoolColumns()函數(shù)每次產(chǎn)生兩個Vec32c類型的向量vWhichCol和 vd,其中vWhichCol對應(yīng)概率矩陣中的目標(biāo)列,算法最終輸出的32個整數(shù)將分別通過掃描這些目標(biāo)列來產(chǎn)生。Sampler()子函數(shù)則依次完成這些掃描任務(wù),且每次掃描時的d值分別對應(yīng)于vd中的分量。在相同實現(xiàn)環(huán)境下的實驗結(jié)果表明,當(dāng)σ=6時,使用上述向量化實現(xiàn)方法可以使得取樣速度提升至14.9×106samples/s。

5 結(jié)束語

本文給出一個整數(shù)上離散高斯取樣的常數(shù)時間實現(xiàn)方法。通過使用單指令多數(shù)據(jù)將其進(jìn)行向量化操作,并在單指令多數(shù)據(jù)的處理器上使用向量化技術(shù),以獲得更好的取樣速度。該方法預(yù)計算過程簡單,存儲量較小,對于存儲資源受限的電子設(shè)備更加友好。本文的實現(xiàn)方法預(yù)計算過程簡單,下一步將研究中心可變的(非中心化)離散高斯分布常數(shù)運行時間的取樣算法,以盡可能減少預(yù)計算存儲量。