Barrett Reduction (1)

Barrett Reduction (1)

はじめに

モジュラー演算で高速に処理したい事案が出てきたので備忘録として残します。

メモ程度なので細かい記載は順次追加していく予定です。

$$a \cdot b \,\, \bmod q$$

$$w = K = \lceil \log_2(q) \rceil$$

$$\mu = \left\lfloor \frac{2^{2K}}{q} \right\rfloor$$

ここで、除算される数のビット長を $L$ とすると、法とする数のビット長 $K$ に対して下式を満たすとき、バレットリダクションの近似計算結果は、厳密に計算したものと一致する。

$$0 \leq L \leq 2K$$

多項式環において、上記が満たされることを今後利用したい。

1. とりあえず実装

下記に注意しながらとりあえず拙く計算例をPythonで実装。

• $w$ は $q$ のビット長になる
• $2^{2w}$ での除算は、シフト演算に相当する

import numpy as np

a = 4
b = 39
q = 40

# normal
print(a*b % q)

# Barret Reduction
w = q.bit_length()
mu = np.floor(2**(2*w) / q).astype(np.int64)

t = a*b
s = np.right_shift(t*mu, 2*w)

c = t - s*q
c = c - q * (c >= q)

2. 修正版 Barrett Reduction ～並列計算に思いを馳せて～

上記をGPU上で計算しようとすると、 $a \cdot b \cdot \mu$ の計算結果が最大 $3K$ ビットになるので、メモリ負荷を軽減したい。

下記のように修正版 Barret Reductionを利用して、最大 $2K$ ビット程度まで改善できる。

import numpy as np

a = 4
b = 39
q = 40

# Barret Reduction Optimized for GPU
w = q.bit_length()
mu = np.floor(2**(2*w) / q).astype(np.int64)

t = a * b
x_1 = np.right_shift(t, w - 2)
x_2 = mu * x_1
s = np.right_shift(x_2, w + 2)
r = s * q
c = t - r
c = c - q * (c >= q)

3. 今後やること

• ビット長が、はるかに大きくなったときへの対応。

参考文献

• https://en.wikipedia.org/wiki/Barrett_reduction