参考文献:
[LS19] 给出了支持 NTT 的 NTRU 变体。他们给出了 AVX2 优化的算法实现,计算效率比 NIST-PQC 所接受的 NTRU 方案要高得多(Gen: 8X, Encap: 5X, Decap: 8X),并且也比其他的 KEM 方案明显更快。
[ACD+18] 使用 [APS15] 的 LWE estimator,评估了提交给 NIST-PQC 的所有 LWE-based、NTRU-based 的 PKE、KEM、Sign 的安全性。对于 128 比特的安全性,环的维度应当介于 700 到 800 之间。然而,radix-2 NTT 仅能处理维度是二的幂次的分圆环,这导致了安全级别的跳跃。NTRU-HRSS 采用了 701 维度的环,NTRU-Prime 采用了 761 维度的环,两者都无法兼容 NTT 算法。而基于 RLWE 的那些方案,可以采取 [BGV12] 的 Generalized LWE 假设将它划分为 3 个维度 256 的小环,从而可以使用 radix-2 NTT 算法。
NTRU-based 相较于 RLWE-based 的优势,
NTRU-based 的缺点:在某些应用中(FHE、可验证加密、群签名)需要模数和噪声规模之间的 gap 足够大。但是 [ABD16] 指出,对于 “过度拉伸”(overstretched)的 NTRU 问题,子域攻击(subfield attack)相当的奏效。并且 NTRU 假设总是比 (R/M)LWE 假设更强(或者采取 [SS11] 的可证明安全变体),并且不同的安全级别下的 NTT 实现无法复用。
所以 NTRU 的参数应当仔细考虑,先要保证方案的安全性,考虑效率才有意义。
我们考虑
m
=
2
k
3
l
,
k
,
l
≥
1
m=2^k3^l,\,\, k,l\ge 1
m=2k3l,k,l≥1 的特殊分圆多项式环,维度
n
=
ϕ
(
m
)
=
2
k
3
l
−
1
=
m
/
3
n=\phi(m)=2^k3^{l-1}=m/3
n=ϕ(m)=2k3l−1=m/3,此时的分圆多项式形如:
Φ
m
(
x
)
=
x
n
−
x
n
/
2
+
1
\Phi_m(x) = x^{n} - x^{n/2} + 1
Φm(x)=xn−xn/2+1
我们考虑
X
2
−
X
+
1
X^2-X+1
X2−X+1 在代数闭包上的分解
(
X
−
α
0
)
⋅
(
X
−
α
1
)
(X-\alpha_0)\cdot(X-\alpha_1)
(X−α0)⋅(X−α1),这要求
α
0
+
α
1
=
1
\alpha_0+\alpha_1=1
α0+α1=1 以及
α
0
α
1
=
1
\alpha_0\alpha_1=1
α0α1=1,我们设置
α
0
=
ζ
6
\alpha_0=\zeta_6
α0=ζ6,
α
1
=
1
−
ζ
6
=
ζ
6
5
\alpha_1=1-\zeta_6=\zeta_6^5
α1=1−ζ6=ζ65,它们都是
6
6
6 次本原单位根。
只要工作域
F
\mathbb F
F 上存在
ζ
6
\zeta_6
ζ6,那么就有如下的环同构:
F
[
x
]
/
(
Φ
m
(
x
)
)
≅
F
[
x
]
/
(
x
n
/
2
−
ζ
6
)
×
F
[
x
]
/
(
x
n
/
2
−
ζ
6
5
)
\mathbb F[x]/(\Phi_m(x)) \cong \mathbb F[x]/(x^{n/2}-\zeta_6) \times\mathbb F[x]/(x^{n/2}-\zeta_6^5)
F[x]/(Φm(x))≅F[x]/(xn/2−ζ6)×F[x]/(xn/2−ζ65)
给定
f
(
x
)
=
∑
i
f
i
x
i
∈
F
[
x
]
/
(
Φ
m
(
x
)
)
f(x)=\sum_i f_i x^i \in \mathbb F[x]/(\Phi_m(x))
f(x)=∑ifixi∈F[x]/(Φm(x)),具体的映射为:
f
0
(
x
)
=
∑
i
=
0
m
/
p
−
1
(
f
i
+
ζ
6
⋅
f
i
+
n
/
2
)
⋅
x
i
f
1
(
x
)
=
∑
i
=
0
m
/
p
−
1
(
f
i
+
f
i
+
n
/
2
−
ζ
6
⋅
f
i
+
n
/
2
)
⋅
x
i
这个蝴蝶的开销是:一次数乘、两次加法、一次减法,基本上和 CT 蝴蝶(一次数乘、一次加法、一次减法)的效率一样。接下来,分别对两个小的卷积环执行 FFT/NTT 即可。
[LS19] 选取的参数:
因此,有以下的环同构:
Z
7681
[
X
]
/
(
X
768
−
X
384
+
1
)
≅
Z
7681
[
X
]
/
(
X
384
−
ζ
6
)
×
Z
7681
[
X
]
/
(
X
384
−
ζ
6
5
)
Z
7681
[
X
]
/
(
X
384
−
ζ
6
)
≅
∏
j
=
0
128
Z
7681
[
X
]
/
(
X
3
−
ζ
768
1
+
j
⋅
6
)
Z
7681
[
X
]
/
(
X
384
−
ζ
6
5
)
≅
∏
j
=
0
128
Z
7681
[
X
]
/
(
X
3
−
ζ
768
5
+
j
⋅
6
)
第一个同构花费
1
1
1 层蝴蝶,其开销接近于 CT 蝴蝶。另外两个同构花费
7
7
7 层蝴蝶,每一层的都是 CT 蝴蝶。共计
8
8
8 层迭代,由于
Z
7681
\mathbb Z_{7681}
Z7681 中不存在
ζ
2304
\zeta_{2304}
ζ2304,所以最终的
deg
=
3
\deg=3
deg=3 的那些 Base Rings 无法继续分解。
InvNTT 可以直接复用 NTT,不过要注意 “比特串翻转” 导致的系数置换。或者采取 GS 蝴蝶,那就不必考虑这些置换,效率可能会更高一些。
定义 modular binomial distribution 简记为 β k \beta_k βk,它的生成过程是:
[LS19] 采取了
β
2
\beta_2
β2 分布,本当
[
1
,
4
,
6
,
4
,
1
]
[1,4,6,4,1]
[1,4,6,4,1],模掉
3
3
3 之后,
Pr
[
−
1
]
=
Pr
[
1
]
=
5
16
,
Pr
[
0
]
=
6
16
\Pr[-1] = \Pr[1] = \dfrac{5}{16},\,\, \Pr[0] = \dfrac{6}{16}
Pr[−1]=Pr[1]=165,Pr[0]=166
如果存在 “随机性恢复”(randomness-recovering)算法
r
←
R
e
c
(
m
,
c
,
p
k
)
r \gets Rec(m,c,pk)
r←Rec(m,c,pk),使得
c
=
E
n
c
(
p
k
,
m
;
r
)
c=Enc(pk,m;r)
c=Enc(pk,m;r),我们称这个 PKE 方案是消息可验证的(message-verifiable)
简单采用 [HPS98] [HS00] 所描述的 NTRU 方案,只不过代数结构从:卷积环 Z [ X ] / ( X N − 1 ) \mathbb Z[X]/(X^N-1) Z[X]/(XN−1),其中的 N N N 是素数(无法使用 FFT/NTT),替换为了:分圆环 Z [ X ] / ( Φ 2304 ( X ) ) \mathbb Z[X]/(\Phi_{2304}(X)) Z[X]/(Φ2304(X)),从而支持 FFT/NTT 算法。

我们额外要求 NTRU 方案的随机性恢复性质:这将使得 IND-CCA2 的归约更加紧致(仅用于此)。由于
c
=
h
r
+
m
c=hr+m
c=hr+m,因此有
R
e
c
(
m
,
c
,
h
)
:
=
(
c
−
m
)
⋅
h
−
1
∈
R
q
Rec(m,c,h) := (c-m) \cdot h^{-1} \in R_q
Rec(m,c,h):=(c−m)⋅h−1∈Rq
这要求公钥
h
=
3
g
/
f
h=3g/f
h=3g/f 也是可逆的。可逆性检查,就是要求 NTT 域的全部系数都是可逆的(非零)。因为
h
h
h 是均匀的,所以
256
256
256 个
Z
q
3
\mathbb Z_q^3
Zq3 上系数当中存在零的概率仅为
Pr
≤
1
−
256
/
q
3
≈
1
−
2
−
30
\Pr \le 1-256/q^3 \approx 1-2^{-30}
Pr≤1−256/q3≈1−2−30,我们完全可以在 PKE 中忽略这个检查,并不对 PKE 方案有实质影响。
OW-CPA PKE 的消息空间、随机带空间:
M
=
R
=
{
f
∈
{
−
1
,
0
,
1
}
768
}
⊆
Z
[
X
]
/
(
Φ
2304
(
X
)
)
M=R=\{f \in \{-1,0,1\}^{768}\} \subseteq \mathbb Z[X]/(\Phi_{2304}(X))
M=R={f∈{−1,0,1}768}⊆Z[X]/(Φ2304(X))
它们的分布记为
D
M
,
D
R
D_M, D_R
DM,DR,解密失败率:
Pr
(
s
k
,
p
k
)
←
G
e
n
,
m
←
D
M
,
r
←
D
R
[
D
e
c
(
s
k
,
E
n
c
(
p
k
,
m
;
r
)
)
≠
m
]
=
ϵ
\underset{(sk,pk)\gets Gen, m\gets D_M, r\gets D_R}{\Pr}[Dec(sk, Enc(pk,m;r)) \neq m] = \epsilon
(sk,pk)←Gen,m←DM,r←DRPr[Dec(sk,Enc(pk,m;r))=m]=ϵ
由于仅当
f
c
=
3
(
g
r
+
f
′
m
)
+
m
fc=3(gr+f'm)+m
fc=3(gr+f′m)+m 的绝对值越过了
q
/
2
q/2
q/2,此时
Z
q
\mathbb Z_q
Zq 无法正确模拟
Z
\mathbb Z
Z,出现解密错误。NIST-POC 接受的 NTRU-HRSS 和 NTRU-Prime,在它们的参数设置下,没有解密错误。[LS19] 分析了
g
,
r
,
f
′
,
m
g,r,f',m
g,r,f′,m 全都服从
β
2
\beta_2
β2 时的解密失败率:在上述 NTTRU 的参数下,错误率为
ϵ
≈
2
−
1230
\epsilon \approx 2^{-1230}
ϵ≈2−1230
但是明文空间的规模 ∣ M ∣ = 3 768 |M|=3^{768} ∣M∣=3768 过于巨大。对于随机选取的 ( s k , p k ) (sk,pk) (sk,pk),存在某消息出现解密失败的概率 Pr ≤ ϵ ⋅ ∣ M ∣ ≈ 2 − 13 \Pr \le \epsilon \cdot|M|\approx 2^{-13} Pr≤ϵ⋅∣M∣≈2−13 相当的大,这可能会导致有效的解密失败攻击(decryption-error attacks)
[LS19] 提出一种转换方法,把消息空间约简到 M ′ = { 0 , 1 } 256 M'=\{0,1\}^{256} M′={0,1}256,使得 ϵ ⋅ ∣ M ∣ ≪ 2 − 128 \epsilon \cdot|M| \ll 2^{-128} ϵ⋅∣M∣≪2−128,从而足够抵御上述攻击,达到 128 比特安全性。

代价是额外的 XOF 和 Hash 的计算,以及密文规模扩大了
256
256
256 比特(对称密文
u
=
E
(
m
′
)
u=E(m')
u=E(m′) 部分)。可以证明,如果存在敌手访问
μ
\mu
μ 次
H
D
M
H_{D_M}
HDM 之后以优势
δ
\delta
δ 打破
C
P
A
′
CPA'
CPA′ 的安全性,那么就存在另一个敌手以优势
δ
−
(
μ
+
1
)
/
∣
M
′
∣
\delta-(\mu+1)/|M'|
δ−(μ+1)/∣M′∣ 打破
C
P
A
CPA
CPA 的安全性,
OW-CPA big
≤
OW-CPA small
\text{ OW-CPA big }\le\text{ OW-CPA small }
OW-CPA big ≤ OW-CPA small
采用 [FO13] [Den02] 的标准提升技术,可以将上述的 OW-CPA PKE 方案,转化为 IND-CCA2 KEM 方案。

归约流程是,
归约链:
OW-CPA usual
≤
OW-CPA message-verifiable
≤
IND-CCA2 KEM
\text{ OW-CPA usual }\le\text{ OW-CPA message-verifiable }\le\text{ IND-CCA2 KEM}
OW-CPA usual ≤ OW-CPA message-verifiable ≤ IND-CCA2 KEM
Hensel remainder:给定模数 q ∈ Z + q \in \mathbb Z^+ q∈Z+ 和字大小 β ∈ Z + \beta \in \mathbb Z^+ β∈Z+,满足 gcd ( q , β ) = 1 \gcd(q,\beta)=1 gcd(q,β)=1 以及 q < β / 2 q<\beta/2 q<β/2。任意的整数 a ∈ Z a \in \mathbb Z a∈Z,可以表示为 a = m q + r β a=mq+r\beta a=mq+rβ,如果限制 m ∈ [ − β / 2 , β / 2 ) m \in [-\beta/2, \beta/2) m∈[−β/2,β/2),那么 r r r 是唯一的。
Montgomery reduction algorithm,记为 M o n t ( a , b ) = a b β − 1 ( m o d q ) Mont(a,b) = ab\beta^{-1}\pmod q Mont(a,b)=abβ−1(modq),
我们选取 β = 2 16 \beta=2^{16} β=216(计算机的半精度),它大于 q = 7681 q=7681 q=7681 的两倍,并且和它互素。同时,step 3 中的关于 β \beta β 的取模、除法,被简化为了 AND、Shift。
AVX2 提供了半精度乘法运算符,给定两个 16 16 16 比特整数 a , b a,b a,b(半精度数),乘积是单精度数(一个字 32 32 32 比特)
因此, M o n t ( a , b ) Mont(a,b) Mont(a,b) 可以优化为:
这导致了更加稠密的 AVX2 向量化(相较于单精度存储),并且节约了一些运算。
现在我们将这个模乘算法,整合到 NTT 算法中,
由于 NTT 蝴蝶中的全部乘法,都是本原根(常数
ζ
\zeta
ζ)和系数(变量
f
f
f)的乘法,因此可以预计算如下的常数:
ζ
′
=
ζ
⋅
β
(
m
o
d
q
)
\zeta'=\zeta \cdot \beta \pmod q
ζ′=ζ⋅β(modq)
那么
M
o
n
t
(
f
,
ζ
′
)
=
f
⋅
ζ
(
m
o
d
q
)
Mont(f,\zeta') = f \cdot \zeta \pmod q
Mont(f,ζ′)=f⋅ζ(modq),这正是我们预期的模乘结果。
[Sei18] 给出的另一个重要的优化是,继续再预计算如下的常数:
ζ
′
′
=
ζ
′
⋅
q
′
(
m
o
d
β
)
\zeta''=\zeta'\cdot q' \pmod \beta
ζ′′=ζ′⋅q′(modβ)
常数特化的
M
o
n
t
c
o
n
s
t
(
f
,
ζ
′
)
Mont_{const}(f,\zeta')
Montconst(f,ζ′) 算法步骤为,
利用预计算的 ζ ′ , ζ ′ ′ \zeta',\zeta'' ζ′,ζ′′,以及 M o n t c o n s t ( f , ζ ′ ) Mont_{const}(f,\zeta') Montconst(f,ζ′),可以计算出正确的 NTT/InvNTT 结果。
对于蝴蝶中的模加运算,可以采取一般性的 Barrett 算法(需要一些乘法)。不过鉴于 q = 7681 q=7681 q=7681 的稀疏比特串表示,[Seil18] 给出了专用的模约减算法(不需要乘法),

对于某一层的蝴蝶(包含若干个多项式,各自分解为长度一半的小多项式),[LS19] 采取了如下的系数打包:
由于 m u l l o , m u l h i mullo, mulhi mullo,mulhi 的计算延迟是 5 5 5-cycles, a d d , s u b add, sub add,sub 的计算延迟是 1 1 1-cycle,从而优化版本的 Montgomery 的延迟是 11 11 11-cycles。为了提高利用率(乱序执行的效果不一定好),[LS19] 手动安排了 6 6 6 个 AVX256 寄存器(加载了 96 96 96 个半精度系数),它们的蝴蝶运算中的乘法是相互独立的,可以填满 CPU 流水线。
为了减少 Load 和 Store 操作(现代处理器的存储墙),[LS19] 采取了层融合技术:长度 768 768 768 的多项式迭代 3 3 3 层分解为 8 8 8 个长度 96 96 96 的多项式。对于每个长度 96 96 96 的多项式(打包在 6 6 6 个 AVX256 寄存器内),持续执行 5 5 5 层 radix-2 NTT(而非只分解一层,store 结果,再 load 下一个多项式),得到 32 32 32 个长度 3 3 3 的最终的分解结果。然后再移动到下一个多项式。
现在我们考虑 Base Ring 上的运算,卷积环 Z q [ X ] / ( X 3 − ζ ) \mathbb Z_q[X]/(X^3-\zeta) Zq[X]/(X3−ζ)
两个元素
f
=
∑
i
f
i
,
g
=
∑
j
g
j
f=\sum_i f_i,\,\, g=\sum_j g_j
f=∑ifi,g=∑jgj 的乘积
h
=
f
g
h=fg
h=fg,写成矩阵形式:
h
=
f
⋅
g
=
f
⋅
g
0
+
f
x
⋅
g
1
+
f
x
2
⋅
g
2
=
[
f
0
ζ
f
2
ζ
f
1
f
1
f
0
ζ
f
2
f
2
f
1
f
0
]
⋅
[
g
0
g
1
g
2
]
=
[
h
0
h
1
h
2
]
使用 Montgomery 算法来计算模乘,
而对于除法 h = g / f h=g/f h=g/f,需要对 Rotation Matrix 求逆:先计算伴随矩阵,然后除以行列式。易知这个逆阵(如果存在的话)也是 Rotation Matrix,
计算伴随矩阵:简单地计算第一列的余子式
f
0
′
=
f
0
2
−
ζ
f
1
f
2
f
1
′
=
ζ
f
2
2
−
f
0
f
1
f
2
′
=
f
1
2
−
f
0
f
2
伴随矩阵
f
∗
=
f
0
′
+
f
1
′
X
+
f
2
′
X
2
f^*=f_0'+f_1'X+f_2'X^2
f∗=f0′+f1′X+f2′X2
计算行列式的逆:简单地计算行列式,可以简化为
d
=
f
0
f
0
′
+
ζ
(
f
1
f
2
′
+
f
2
f
1
′
)
d = f_0f_0' + \zeta(f_1f_2' + f_2f_1')
d=f0f0′+ζ(f1f2′+f2f1′)
然后计算
d
−
1
=
d
q
−
2
d^{-1}=d^{q-2}
d−1=dq−2,利用快速幂算法
最后得到 f − 1 = d − 1 f ∗ f^{-1}=d^{-1}f^* f−1=d−1f∗,其中的所有模乘运算都采取 M o n t ( f , g ) Mont(f,g) Mont(f,g) 和 M o n t c o n s t ( f , ζ ′ ) Mont_{const}(f,\zeta') Montconst(f,ζ′) 来计算
由于 M o n t ( f , g ) Mont(f,g) Mont(f,g) 会引入一些因子 ( β − 1 ) v (\beta^{-1})^v (β−1)v,我们需要追踪它们,结果是:
也使用 AVX2 实现这些 Base Ring 上的运算。那么环 Z 7681 [ X ] / ( X 768 − X 384 + 1 ) \mathbb Z_{7681}[X]/(X^{768}-X^{384}+1) Z7681[X]/(X768−X384+1) 上的多项式运算:NTT,InvNTT,小环 Z 7681 [ X ] / ( X 3 − ζ ) \mathbb Z_{7681}[X]/(X^3-\zeta) Z7681[X]/(X3−ζ) 上的乘法、除法,它们的计算速度都是极快的。这导致效率瓶颈,反而是 Sample、Hash、XOF、Pack/Unpack 这些运算。
为了快速采样 β 2 = ( a 1 + a 2 ) − ( b 1 + b 2 ) ( m o d ± 3 ) \beta_2=(a_1+a_2)-(b_1+b_2)\pmod{^\pm3} β2=(a1+a2)−(b1+b2)(mod±3),[LS19] 采取查表的方法。构造如下的 LUT,它是长度为 16 16 16 的向量:
| ( a 1 a 2 b 1 b 2 ) 2 (a_1a_2b_1b_2)_2 (a1a2b1b2)2 | β 2 \beta_2 β2 |
|---|---|
| 0 = 0000 0=0000 0=0000 | 0 0 0 |
| 1 = 0001 1=0001 1=0001 | − 1 -1 −1 |
| 2 = 0010 2=0010 2=0010 | − 1 -1 −1 |
| 3 = 0011 3=0011 3=0011 | 1 1 1 |
| … | … |
| 14 = 1110 14=1110 14=1110 | 1 1 1 |
| 15 = 1111 15=1111 15=1111 | 0 0 0 |
由于
β
2
\beta_2
β2 取值
{
−
1
,
0
,
1
}
\{-1,0,1\}
{−1,0,1} 只花费
2
2
2 比特,因此上述的 LUT 可以被存储为一个字
T
T
T,采样算法就是:产生均匀的四比特
i
=
(
a
1
a
2
b
1
b
2
)
2
i=(a_1a_2b_1b_2)_2
i=(a1a2b1b2)2,然后简单的 shift,
β
2
[
a
1
,
a
2
,
b
2
,
b
2
]
=
(
T
≫
(
2
i
)
)
[
0
]
\beta_2[a_1,a_2,b_2,b_2] = \left(T \gg (2i)\right)[0]
β2[a1,a2,b2,b2]=(T≫(2i))[0]
考虑到表格的对称性
β
2
[
i
]
=
−
β
2
[
15
−
i
]
\beta_2[i]=-\beta_2[15-i]
β2[i]=−β2[15−i],实际上 LUT 可以被存储为半个字(
16
16
16 比特),从而一个 AVX256 寄存器中可以填入
16
16
16 张表,以
16
16
16 路并行的方式快速采样。不过 AVX2 不支持半精度整数的移位,所以还需要把 LUT 重排一下。
为了生成 f , g , r f,g,r f,g,r,我们使用 XOF 对某个 seed 做扩展,然后调用上述的 Binomial Sampler 执行各个系数的采样。但是 SHAKE 的速度很慢(相对于 NTT 而言),所以 [LS19] 选用了 CTR-mode AES 作为 XOF
在 KEM 的封装算法中,需要计算两个哈希 r ← H D R ( m ) r \gets H_{D_R}(m) r←HDR(m) 和 k ← H K ( m ) k \gets H_{K}(m) k←HK(m),[LS19] 使用 SHA512 合并地计算出 512 512 512 比特的摘要,然后各自设置 r , k ∈ { 0 , 1 } 256 r,k \in \{0,1\}^{256} r,k∈{0,1}256 是它的各一半。
由于 q = 7681 ≈ 2 13 q=7681 \approx 2^{13} q=7681≈213,因此我们将多项式的每 8 8 8 个系数(本来占据 16 16 16 字节),打包到连续的 13 13 13 字节,减少通信开销。这个打包过程也是极慢的,如果不进行 AVX2 优化,它甚至比 SHA3 的速度还要慢。
鉴于 AVX256 寄存器是连续读取,[LS19] 间隔 16 16 16 提取 8 8 8 个系数(而非连续的 8 8 8 个系数),将它们打包在一起。于是可以将 16 × 8 16 \times 8 16×8 个系数,连续读取到 8 8 8 个 AVX256 寄存器中(每个寄存器加载连续的 16 16 16 个系数),然后并行地打包,得到 16 16 16 个长度 13 13 13 字节的打包结果。