1*0957b409SSimon J. Gerraty /*
2*0957b409SSimon J. Gerraty * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org>
3*0957b409SSimon J. Gerraty *
4*0957b409SSimon J. Gerraty * Permission is hereby granted, free of charge, to any person obtaining
5*0957b409SSimon J. Gerraty * a copy of this software and associated documentation files (the
6*0957b409SSimon J. Gerraty * "Software"), to deal in the Software without restriction, including
7*0957b409SSimon J. Gerraty * without limitation the rights to use, copy, modify, merge, publish,
8*0957b409SSimon J. Gerraty * distribute, sublicense, and/or sell copies of the Software, and to
9*0957b409SSimon J. Gerraty * permit persons to whom the Software is furnished to do so, subject to
10*0957b409SSimon J. Gerraty * the following conditions:
11*0957b409SSimon J. Gerraty *
12*0957b409SSimon J. Gerraty * The above copyright notice and this permission notice shall be
13*0957b409SSimon J. Gerraty * included in all copies or substantial portions of the Software.
14*0957b409SSimon J. Gerraty *
15*0957b409SSimon J. Gerraty * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
16*0957b409SSimon J. Gerraty * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
17*0957b409SSimon J. Gerraty * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
18*0957b409SSimon J. Gerraty * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS
19*0957b409SSimon J. Gerraty * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN
20*0957b409SSimon J. Gerraty * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN
21*0957b409SSimon J. Gerraty * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
22*0957b409SSimon J. Gerraty * SOFTWARE.
23*0957b409SSimon J. Gerraty */
24*0957b409SSimon J. Gerraty
25*0957b409SSimon J. Gerraty #include "inner.h"
26*0957b409SSimon J. Gerraty
27*0957b409SSimon J. Gerraty /*
28*0957b409SSimon J. Gerraty * This implementation uses 32-bit multiplications, and only the low
29*0957b409SSimon J. Gerraty * 32 bits for each multiplication result. This is meant primarily for
30*0957b409SSimon J. Gerraty * the ARM Cortex M0 and M0+, whose multiplication opcode does not yield
31*0957b409SSimon J. Gerraty * the upper 32 bits; but it might also be useful on architectures where
32*0957b409SSimon J. Gerraty * access to the upper 32 bits requires use of specific registers that
33*0957b409SSimon J. Gerraty * create contention (e.g. on i386, "mul" necessarily outputs the result
34*0957b409SSimon J. Gerraty * in edx:eax, while "imul" can use any registers but is limited to the
35*0957b409SSimon J. Gerraty * low 32 bits).
36*0957b409SSimon J. Gerraty *
37*0957b409SSimon J. Gerraty * The implementation trick that is used here is bit-reversing (bit 0
38*0957b409SSimon J. Gerraty * is swapped with bit 31, bit 1 with bit 30, and so on). In GF(2)[X],
39*0957b409SSimon J. Gerraty * for all values x and y, we have:
40*0957b409SSimon J. Gerraty * rev32(x) * rev32(y) = rev64(x * y)
41*0957b409SSimon J. Gerraty * In other words, if we bit-reverse (over 32 bits) the operands, then we
42*0957b409SSimon J. Gerraty * bit-reverse (over 64 bits) the result.
43*0957b409SSimon J. Gerraty */
44*0957b409SSimon J. Gerraty
45*0957b409SSimon J. Gerraty /*
46*0957b409SSimon J. Gerraty * Multiplication in GF(2)[X], truncated to its low 32 bits.
47*0957b409SSimon J. Gerraty */
48*0957b409SSimon J. Gerraty static inline uint32_t
bmul32(uint32_t x,uint32_t y)49*0957b409SSimon J. Gerraty bmul32(uint32_t x, uint32_t y)
50*0957b409SSimon J. Gerraty {
51*0957b409SSimon J. Gerraty uint32_t x0, x1, x2, x3;
52*0957b409SSimon J. Gerraty uint32_t y0, y1, y2, y3;
53*0957b409SSimon J. Gerraty uint32_t z0, z1, z2, z3;
54*0957b409SSimon J. Gerraty
55*0957b409SSimon J. Gerraty x0 = x & (uint32_t)0x11111111;
56*0957b409SSimon J. Gerraty x1 = x & (uint32_t)0x22222222;
57*0957b409SSimon J. Gerraty x2 = x & (uint32_t)0x44444444;
58*0957b409SSimon J. Gerraty x3 = x & (uint32_t)0x88888888;
59*0957b409SSimon J. Gerraty y0 = y & (uint32_t)0x11111111;
60*0957b409SSimon J. Gerraty y1 = y & (uint32_t)0x22222222;
61*0957b409SSimon J. Gerraty y2 = y & (uint32_t)0x44444444;
62*0957b409SSimon J. Gerraty y3 = y & (uint32_t)0x88888888;
63*0957b409SSimon J. Gerraty z0 = (x0 * y0) ^ (x1 * y3) ^ (x2 * y2) ^ (x3 * y1);
64*0957b409SSimon J. Gerraty z1 = (x0 * y1) ^ (x1 * y0) ^ (x2 * y3) ^ (x3 * y2);
65*0957b409SSimon J. Gerraty z2 = (x0 * y2) ^ (x1 * y1) ^ (x2 * y0) ^ (x3 * y3);
66*0957b409SSimon J. Gerraty z3 = (x0 * y3) ^ (x1 * y2) ^ (x2 * y1) ^ (x3 * y0);
67*0957b409SSimon J. Gerraty z0 &= (uint32_t)0x11111111;
68*0957b409SSimon J. Gerraty z1 &= (uint32_t)0x22222222;
69*0957b409SSimon J. Gerraty z2 &= (uint32_t)0x44444444;
70*0957b409SSimon J. Gerraty z3 &= (uint32_t)0x88888888;
71*0957b409SSimon J. Gerraty return z0 | z1 | z2 | z3;
72*0957b409SSimon J. Gerraty }
73*0957b409SSimon J. Gerraty
74*0957b409SSimon J. Gerraty /*
75*0957b409SSimon J. Gerraty * Bit-reverse a 32-bit word.
76*0957b409SSimon J. Gerraty */
77*0957b409SSimon J. Gerraty static uint32_t
rev32(uint32_t x)78*0957b409SSimon J. Gerraty rev32(uint32_t x)
79*0957b409SSimon J. Gerraty {
80*0957b409SSimon J. Gerraty #define RMS(m, s) do { \
81*0957b409SSimon J. Gerraty x = ((x & (uint32_t)(m)) << (s)) \
82*0957b409SSimon J. Gerraty | ((x >> (s)) & (uint32_t)(m)); \
83*0957b409SSimon J. Gerraty } while (0)
84*0957b409SSimon J. Gerraty
85*0957b409SSimon J. Gerraty RMS(0x55555555, 1);
86*0957b409SSimon J. Gerraty RMS(0x33333333, 2);
87*0957b409SSimon J. Gerraty RMS(0x0F0F0F0F, 4);
88*0957b409SSimon J. Gerraty RMS(0x00FF00FF, 8);
89*0957b409SSimon J. Gerraty return (x << 16) | (x >> 16);
90*0957b409SSimon J. Gerraty
91*0957b409SSimon J. Gerraty #undef RMS
92*0957b409SSimon J. Gerraty }
93*0957b409SSimon J. Gerraty
94*0957b409SSimon J. Gerraty /* see bearssl_hash.h */
95*0957b409SSimon J. Gerraty void
br_ghash_ctmul32(void * y,const void * h,const void * data,size_t len)96*0957b409SSimon J. Gerraty br_ghash_ctmul32(void *y, const void *h, const void *data, size_t len)
97*0957b409SSimon J. Gerraty {
98*0957b409SSimon J. Gerraty /*
99*0957b409SSimon J. Gerraty * This implementation is similar to br_ghash_ctmul() except
100*0957b409SSimon J. Gerraty * that we have to do the multiplication twice, with the
101*0957b409SSimon J. Gerraty * "normal" and "bit reversed" operands. Hence we end up with
102*0957b409SSimon J. Gerraty * eighteen 32-bit multiplications instead of nine.
103*0957b409SSimon J. Gerraty */
104*0957b409SSimon J. Gerraty
105*0957b409SSimon J. Gerraty const unsigned char *buf, *hb;
106*0957b409SSimon J. Gerraty unsigned char *yb;
107*0957b409SSimon J. Gerraty uint32_t yw[4];
108*0957b409SSimon J. Gerraty uint32_t hw[4], hwr[4];
109*0957b409SSimon J. Gerraty
110*0957b409SSimon J. Gerraty buf = data;
111*0957b409SSimon J. Gerraty yb = y;
112*0957b409SSimon J. Gerraty hb = h;
113*0957b409SSimon J. Gerraty yw[3] = br_dec32be(yb);
114*0957b409SSimon J. Gerraty yw[2] = br_dec32be(yb + 4);
115*0957b409SSimon J. Gerraty yw[1] = br_dec32be(yb + 8);
116*0957b409SSimon J. Gerraty yw[0] = br_dec32be(yb + 12);
117*0957b409SSimon J. Gerraty hw[3] = br_dec32be(hb);
118*0957b409SSimon J. Gerraty hw[2] = br_dec32be(hb + 4);
119*0957b409SSimon J. Gerraty hw[1] = br_dec32be(hb + 8);
120*0957b409SSimon J. Gerraty hw[0] = br_dec32be(hb + 12);
121*0957b409SSimon J. Gerraty hwr[3] = rev32(hw[3]);
122*0957b409SSimon J. Gerraty hwr[2] = rev32(hw[2]);
123*0957b409SSimon J. Gerraty hwr[1] = rev32(hw[1]);
124*0957b409SSimon J. Gerraty hwr[0] = rev32(hw[0]);
125*0957b409SSimon J. Gerraty while (len > 0) {
126*0957b409SSimon J. Gerraty const unsigned char *src;
127*0957b409SSimon J. Gerraty unsigned char tmp[16];
128*0957b409SSimon J. Gerraty int i;
129*0957b409SSimon J. Gerraty uint32_t a[18], b[18], c[18];
130*0957b409SSimon J. Gerraty uint32_t d0, d1, d2, d3, d4, d5, d6, d7;
131*0957b409SSimon J. Gerraty uint32_t zw[8];
132*0957b409SSimon J. Gerraty
133*0957b409SSimon J. Gerraty if (len >= 16) {
134*0957b409SSimon J. Gerraty src = buf;
135*0957b409SSimon J. Gerraty buf += 16;
136*0957b409SSimon J. Gerraty len -= 16;
137*0957b409SSimon J. Gerraty } else {
138*0957b409SSimon J. Gerraty memcpy(tmp, buf, len);
139*0957b409SSimon J. Gerraty memset(tmp + len, 0, (sizeof tmp) - len);
140*0957b409SSimon J. Gerraty src = tmp;
141*0957b409SSimon J. Gerraty len = 0;
142*0957b409SSimon J. Gerraty }
143*0957b409SSimon J. Gerraty yw[3] ^= br_dec32be(src);
144*0957b409SSimon J. Gerraty yw[2] ^= br_dec32be(src + 4);
145*0957b409SSimon J. Gerraty yw[1] ^= br_dec32be(src + 8);
146*0957b409SSimon J. Gerraty yw[0] ^= br_dec32be(src + 12);
147*0957b409SSimon J. Gerraty
148*0957b409SSimon J. Gerraty /*
149*0957b409SSimon J. Gerraty * We are using Karatsuba: the 128x128 multiplication is
150*0957b409SSimon J. Gerraty * reduced to three 64x64 multiplications, hence nine
151*0957b409SSimon J. Gerraty * 32x32 multiplications. With the bit-reversal trick,
152*0957b409SSimon J. Gerraty * we have to perform 18 32x32 multiplications.
153*0957b409SSimon J. Gerraty */
154*0957b409SSimon J. Gerraty
155*0957b409SSimon J. Gerraty /*
156*0957b409SSimon J. Gerraty * y[0,1]*h[0,1] -> 0,1,4
157*0957b409SSimon J. Gerraty * y[2,3]*h[2,3] -> 2,3,5
158*0957b409SSimon J. Gerraty * (y[0,1]+y[2,3])*(h[0,1]+h[2,3]) -> 6,7,8
159*0957b409SSimon J. Gerraty */
160*0957b409SSimon J. Gerraty
161*0957b409SSimon J. Gerraty a[0] = yw[0];
162*0957b409SSimon J. Gerraty a[1] = yw[1];
163*0957b409SSimon J. Gerraty a[2] = yw[2];
164*0957b409SSimon J. Gerraty a[3] = yw[3];
165*0957b409SSimon J. Gerraty a[4] = a[0] ^ a[1];
166*0957b409SSimon J. Gerraty a[5] = a[2] ^ a[3];
167*0957b409SSimon J. Gerraty a[6] = a[0] ^ a[2];
168*0957b409SSimon J. Gerraty a[7] = a[1] ^ a[3];
169*0957b409SSimon J. Gerraty a[8] = a[6] ^ a[7];
170*0957b409SSimon J. Gerraty
171*0957b409SSimon J. Gerraty a[ 9] = rev32(yw[0]);
172*0957b409SSimon J. Gerraty a[10] = rev32(yw[1]);
173*0957b409SSimon J. Gerraty a[11] = rev32(yw[2]);
174*0957b409SSimon J. Gerraty a[12] = rev32(yw[3]);
175*0957b409SSimon J. Gerraty a[13] = a[ 9] ^ a[10];
176*0957b409SSimon J. Gerraty a[14] = a[11] ^ a[12];
177*0957b409SSimon J. Gerraty a[15] = a[ 9] ^ a[11];
178*0957b409SSimon J. Gerraty a[16] = a[10] ^ a[12];
179*0957b409SSimon J. Gerraty a[17] = a[15] ^ a[16];
180*0957b409SSimon J. Gerraty
181*0957b409SSimon J. Gerraty b[0] = hw[0];
182*0957b409SSimon J. Gerraty b[1] = hw[1];
183*0957b409SSimon J. Gerraty b[2] = hw[2];
184*0957b409SSimon J. Gerraty b[3] = hw[3];
185*0957b409SSimon J. Gerraty b[4] = b[0] ^ b[1];
186*0957b409SSimon J. Gerraty b[5] = b[2] ^ b[3];
187*0957b409SSimon J. Gerraty b[6] = b[0] ^ b[2];
188*0957b409SSimon J. Gerraty b[7] = b[1] ^ b[3];
189*0957b409SSimon J. Gerraty b[8] = b[6] ^ b[7];
190*0957b409SSimon J. Gerraty
191*0957b409SSimon J. Gerraty b[ 9] = hwr[0];
192*0957b409SSimon J. Gerraty b[10] = hwr[1];
193*0957b409SSimon J. Gerraty b[11] = hwr[2];
194*0957b409SSimon J. Gerraty b[12] = hwr[3];
195*0957b409SSimon J. Gerraty b[13] = b[ 9] ^ b[10];
196*0957b409SSimon J. Gerraty b[14] = b[11] ^ b[12];
197*0957b409SSimon J. Gerraty b[15] = b[ 9] ^ b[11];
198*0957b409SSimon J. Gerraty b[16] = b[10] ^ b[12];
199*0957b409SSimon J. Gerraty b[17] = b[15] ^ b[16];
200*0957b409SSimon J. Gerraty
201*0957b409SSimon J. Gerraty for (i = 0; i < 18; i ++) {
202*0957b409SSimon J. Gerraty c[i] = bmul32(a[i], b[i]);
203*0957b409SSimon J. Gerraty }
204*0957b409SSimon J. Gerraty
205*0957b409SSimon J. Gerraty c[4] ^= c[0] ^ c[1];
206*0957b409SSimon J. Gerraty c[5] ^= c[2] ^ c[3];
207*0957b409SSimon J. Gerraty c[8] ^= c[6] ^ c[7];
208*0957b409SSimon J. Gerraty
209*0957b409SSimon J. Gerraty c[13] ^= c[ 9] ^ c[10];
210*0957b409SSimon J. Gerraty c[14] ^= c[11] ^ c[12];
211*0957b409SSimon J. Gerraty c[17] ^= c[15] ^ c[16];
212*0957b409SSimon J. Gerraty
213*0957b409SSimon J. Gerraty /*
214*0957b409SSimon J. Gerraty * y[0,1]*h[0,1] -> 0,9^4,1^13,10
215*0957b409SSimon J. Gerraty * y[2,3]*h[2,3] -> 2,11^5,3^14,12
216*0957b409SSimon J. Gerraty * (y[0,1]+y[2,3])*(h[0,1]+h[2,3]) -> 6,15^8,7^17,16
217*0957b409SSimon J. Gerraty */
218*0957b409SSimon J. Gerraty d0 = c[0];
219*0957b409SSimon J. Gerraty d1 = c[4] ^ (rev32(c[9]) >> 1);
220*0957b409SSimon J. Gerraty d2 = c[1] ^ c[0] ^ c[2] ^ c[6] ^ (rev32(c[13]) >> 1);
221*0957b409SSimon J. Gerraty d3 = c[4] ^ c[5] ^ c[8]
222*0957b409SSimon J. Gerraty ^ (rev32(c[10] ^ c[9] ^ c[11] ^ c[15]) >> 1);
223*0957b409SSimon J. Gerraty d4 = c[2] ^ c[1] ^ c[3] ^ c[7]
224*0957b409SSimon J. Gerraty ^ (rev32(c[13] ^ c[14] ^ c[17]) >> 1);
225*0957b409SSimon J. Gerraty d5 = c[5] ^ (rev32(c[11] ^ c[10] ^ c[12] ^ c[16]) >> 1);
226*0957b409SSimon J. Gerraty d6 = c[3] ^ (rev32(c[14]) >> 1);
227*0957b409SSimon J. Gerraty d7 = rev32(c[12]) >> 1;
228*0957b409SSimon J. Gerraty
229*0957b409SSimon J. Gerraty zw[0] = d0 << 1;
230*0957b409SSimon J. Gerraty zw[1] = (d1 << 1) | (d0 >> 31);
231*0957b409SSimon J. Gerraty zw[2] = (d2 << 1) | (d1 >> 31);
232*0957b409SSimon J. Gerraty zw[3] = (d3 << 1) | (d2 >> 31);
233*0957b409SSimon J. Gerraty zw[4] = (d4 << 1) | (d3 >> 31);
234*0957b409SSimon J. Gerraty zw[5] = (d5 << 1) | (d4 >> 31);
235*0957b409SSimon J. Gerraty zw[6] = (d6 << 1) | (d5 >> 31);
236*0957b409SSimon J. Gerraty zw[7] = (d7 << 1) | (d6 >> 31);
237*0957b409SSimon J. Gerraty
238*0957b409SSimon J. Gerraty for (i = 0; i < 4; i ++) {
239*0957b409SSimon J. Gerraty uint32_t lw;
240*0957b409SSimon J. Gerraty
241*0957b409SSimon J. Gerraty lw = zw[i];
242*0957b409SSimon J. Gerraty zw[i + 4] ^= lw ^ (lw >> 1) ^ (lw >> 2) ^ (lw >> 7);
243*0957b409SSimon J. Gerraty zw[i + 3] ^= (lw << 31) ^ (lw << 30) ^ (lw << 25);
244*0957b409SSimon J. Gerraty }
245*0957b409SSimon J. Gerraty memcpy(yw, zw + 4, sizeof yw);
246*0957b409SSimon J. Gerraty }
247*0957b409SSimon J. Gerraty br_enc32be(yb, yw[3]);
248*0957b409SSimon J. Gerraty br_enc32be(yb + 4, yw[2]);
249*0957b409SSimon J. Gerraty br_enc32be(yb + 8, yw[1]);
250*0957b409SSimon J. Gerraty br_enc32be(yb + 12, yw[0]);
251*0957b409SSimon J. Gerraty }
252