ntt.mx raw
1 package hamcrypto
2
3 // Number-Theoretic Transform (NTT) over Z_257 for Hamadryad.
4 // n=64, p=257. psi=9 (primitive 128th root of unity).
5
6 // mod257 reduces x into [0, 256].
7 func mod257(x int32) (r uint16) {
8 x %= HamP
9 if x < 0 {
10 x += HamP
11 }
12 return uint16(x)
13 }
14
15 // powMod computes base^exp mod 257.
16 func powMod(base, exp int32) (r uint16) {
17 result := int32(1)
18 b := base % HamP
19 if b < 0 {
20 b += HamP
21 }
22 e := exp
23 for e > 0 {
24 if e&1 == 1 {
25 result = (result * b) % HamP
26 }
27 b = (b * b) % HamP
28 e >>= 1
29 }
30 return uint16(result)
31 }
32
33 // invMod computes the modular inverse of a mod 257.
34 func invMod(a uint16) (r uint16) {
35 return powMod(int32(a), HamP-2)
36 }
37
38 // bitRev6 reverses the low 6 bits of x.
39 func bitRev6(x int32) (r int32) {
40 r = 0
41 xv := x
42 for i := int32(0); i < 6; i++ {
43 r = (r << 1) | (xv & 1)
44 xv >>= 1
45 }
46 return r
47 }
48
49 // Pre-computed tables, filled by initNTTTables.
50 var psiPows [128]uint16
51 var psiInvPows [128]uint16
52 var invN uint16
53 var nttTablesReady bool
54
55 // initNTTTables populates the NTT lookup tables.
56 func initNTTTables() {
57 if nttTablesReady {
58 return
59 }
60 const psi = 9
61
62 psiPows[0] = 1
63 for i := int32(1); i < 128; i++ {
64 psiPows[i] = mod257(int32(psiPows[i-1]) * psi)
65 }
66
67 psiInv := invMod(psi)
68 psiInvPows[0] = 1
69 for i := int32(1); i < 128; i++ {
70 psiInvPows[i] = mod257(int32(psiInvPows[i-1]) * int32(psiInv))
71 }
72
73 invN = invMod(uint16(HamN))
74 nttTablesReady = true
75 }
76
77 // ntt64 computes the forward negacyclic NTT of a length-64 polynomial over Z_257.
78 func ntt64(a *[HamN]uint16) {
79 const n = int32(HamN)
80
81 // Pre-multiply by psi^i.
82 for i := int32(0); i < n; i++ {
83 a[i] = mod257(int32(a[i]) * int32(psiPows[i]))
84 }
85
86 // Bit-reversal permutation.
87 for i := int32(0); i < n; i++ {
88 j := bitRev6(i)
89 if i < j {
90 a[i], a[j] = a[j], a[i]
91 }
92 }
93
94 // Cooley-Tukey butterfly stages.
95 length := int32(1)
96 for length < n {
97 step := n / (2 * length)
98 start := int32(0)
99 for start < n {
100 for j := int32(0); j < length; j++ {
101 tw := psiPows[(2*j*step)%128]
102 u := a[start+j]
103 v := mod257(int32(a[start+j+length]) * int32(tw))
104 a[start+j] = mod257(int32(u) + int32(v))
105 a[start+j+length] = mod257(int32(u) - int32(v))
106 }
107 start += 2 * length
108 }
109 length <<= 1
110 }
111 }
112
113 // intt64 computes the inverse negacyclic NTT, recovering coefficients.
114 func intt64(a *[HamN]uint16) {
115 const n = int32(HamN)
116
117 // Bit-reversal permutation.
118 for i := int32(0); i < n; i++ {
119 j := bitRev6(i)
120 if i < j {
121 a[i], a[j] = a[j], a[i]
122 }
123 }
124
125 // Gentleman-Sande butterfly (inverse DIT).
126 length := int32(1)
127 for length < n {
128 step := n / (2 * length)
129 start := int32(0)
130 for start < n {
131 for j := int32(0); j < length; j++ {
132 tw := psiInvPows[(2*j*step)%128]
133 u := a[start+j]
134 v := mod257(int32(a[start+j+length]) * int32(tw))
135 a[start+j] = mod257(int32(u) + int32(v))
136 a[start+j+length] = mod257(int32(u) - int32(v))
137 }
138 start += 2 * length
139 }
140 length <<= 1
141 }
142
143 // Multiply by 1/n.
144 for i := int32(0); i < n; i++ {
145 a[i] = mod257(int32(a[i]) * int32(invN))
146 }
147
148 // Undo psi pre-multiplication.
149 for i := int32(0); i < n; i++ {
150 a[i] = mod257(int32(a[i]) * int32(psiInvPows[i]))
151 }
152 }
153