mozilla
/
libprio
зеркало из https://github.com/mozilla/libprio.git
1
0
Форкнуть 0
Код Задачи Пакеты Проекты Релизы Вики Активность
libprio/mpi/mp_gf2m.c

679 строки
16 KiB
C
Исходник Постоянная ссылка Ответственный История

/* This Source Code Form is subject to the terms of the Mozilla Public
* License, v. 2.0. If a copy of the MPL was not distributed with this
* file, You can obtain one at http://mozilla.org/MPL/2.0/. */
#include "mp_gf2m.h"
#include "mp_gf2m-priv.h"
#include "mplogic.h"
#include "mpi-priv.h"
const mp_digit mp_gf2m_sqr_tb[16] =
{
0, 1, 4, 5, 16, 17, 20, 21,
64, 65, 68, 69, 80, 81, 84, 85
};
/* Multiply two binary polynomials mp_digits a, b.
* Result is a polynomial with degree < 2 * MP_DIGIT_BITS - 1.
* Output in two mp_digits rh, rl.
*/
#if MP_DIGIT_BITS == 32
void
s_bmul_1x1(mp_digit *rh, mp_digit *rl, const mp_digit a, const mp_digit b)
{
register mp_digit h, l, s;
mp_digit tab[8], top2b = a >> 30;
register mp_digit a1, a2, a4;
a1 = a & (0x3FFFFFFF);
a2 = a1 << 1;
a4 = a2 << 1;
tab[0] = 0;
tab[1] = a1;
tab[2] = a2;
tab[3] = a1 ^ a2;
tab[4] = a4;
tab[5] = a1 ^ a4;
tab[6] = a2 ^ a4;
tab[7] = a1 ^ a2 ^ a4;
s = tab[b & 0x7];
l = s;
s = tab[b >> 3 & 0x7];
l ^= s << 3;
h = s >> 29;
s = tab[b >> 6 & 0x7];
l ^= s << 6;
h ^= s >> 26;
s = tab[b >> 9 & 0x7];
l ^= s << 9;
h ^= s >> 23;
s = tab[b >> 12 & 0x7];
l ^= s << 12;
h ^= s >> 20;
s = tab[b >> 15 & 0x7];
l ^= s << 15;
h ^= s >> 17;
s = tab[b >> 18 & 0x7];
l ^= s << 18;
h ^= s >> 14;
s = tab[b >> 21 & 0x7];
l ^= s << 21;
h ^= s >> 11;
s = tab[b >> 24 & 0x7];
l ^= s << 24;
h ^= s >> 8;
s = tab[b >> 27 & 0x7];
l ^= s << 27;
h ^= s >> 5;
s = tab[b >> 30];
l ^= s << 30;
h ^= s >> 2;
/* compensate for the top two bits of a */
if (top2b & 01) {
l ^= b << 30;
h ^= b >> 2;
}
if (top2b & 02) {
l ^= b << 31;
h ^= b >> 1;
}
*rh = h;
*rl = l;
}
#else
void
s_bmul_1x1(mp_digit *rh, mp_digit *rl, const mp_digit a, const mp_digit b)
{
register mp_digit h, l, s;
mp_digit tab[16], top3b = a >> 61;
register mp_digit a1, a2, a4, a8;
a1 = a & (0x1FFFFFFFFFFFFFFFULL);
a2 = a1 << 1;
a4 = a2 << 1;
a8 = a4 << 1;
tab[0] = 0;
tab[1] = a1;
tab[2] = a2;
tab[3] = a1 ^ a2;
tab[4] = a4;
tab[5] = a1 ^ a4;
tab[6] = a2 ^ a4;
tab[7] = a1 ^ a2 ^ a4;
tab[8] = a8;
tab[9] = a1 ^ a8;
tab[10] = a2 ^ a8;
tab[11] = a1 ^ a2 ^ a8;
tab[12] = a4 ^ a8;
tab[13] = a1 ^ a4 ^ a8;
tab[14] = a2 ^ a4 ^ a8;
tab[15] = a1 ^ a2 ^ a4 ^ a8;
s = tab[b & 0xF];
l = s;
s = tab[b >> 4 & 0xF];
l ^= s << 4;
h = s >> 60;
s = tab[b >> 8 & 0xF];
l ^= s << 8;
h ^= s >> 56;
s = tab[b >> 12 & 0xF];
l ^= s << 12;
h ^= s >> 52;
s = tab[b >> 16 & 0xF];
l ^= s << 16;
h ^= s >> 48;
s = tab[b >> 20 & 0xF];
l ^= s << 20;
h ^= s >> 44;
s = tab[b >> 24 & 0xF];
l ^= s << 24;
h ^= s >> 40;
s = tab[b >> 28 & 0xF];
l ^= s << 28;
h ^= s >> 36;
s = tab[b >> 32 & 0xF];
l ^= s << 32;
h ^= s >> 32;
s = tab[b >> 36 & 0xF];
l ^= s << 36;
h ^= s >> 28;
s = tab[b >> 40 & 0xF];
l ^= s << 40;
h ^= s >> 24;
s = tab[b >> 44 & 0xF];
l ^= s << 44;
h ^= s >> 20;
s = tab[b >> 48 & 0xF];
l ^= s << 48;
h ^= s >> 16;
s = tab[b >> 52 & 0xF];
l ^= s << 52;
h ^= s >> 12;
s = tab[b >> 56 & 0xF];
l ^= s << 56;
h ^= s >> 8;
s = tab[b >> 60];
l ^= s << 60;
h ^= s >> 4;
/* compensate for the top three bits of a */
if (top3b & 01) {
l ^= b << 61;
h ^= b >> 3;
}
if (top3b & 02) {
l ^= b << 62;
h ^= b >> 2;
}
if (top3b & 04) {
l ^= b << 63;
h ^= b >> 1;
}
*rh = h;
*rl = l;
}
#endif
/* Compute xor-multiply of two binary polynomials (a1, a0) x (b1, b0)
* result is a binary polynomial in 4 mp_digits r[4].
* The caller MUST ensure that r has the right amount of space allocated.
*/
void
s_bmul_2x2(mp_digit *r, const mp_digit a1, const mp_digit a0, const mp_digit b1,
const mp_digit b0)
{
mp_digit m1, m0;
/* r[3] = h1, r[2] = h0; r[1] = l1; r[0] = l0 */
s_bmul_1x1(r + 3, r + 2, a1, b1);
s_bmul_1x1(r + 1, r, a0, b0);
s_bmul_1x1(&m1, &m0, a0 ^ a1, b0 ^ b1);
/* Correction on m1 ^= l1 ^ h1; m0 ^= l0 ^ h0; */
r[2] ^= m1 ^ r[1] ^ r[3]; /* h0 ^= m1 ^ l1 ^ h1; */
r[1] = r[3] ^ r[2] ^ r[0] ^ m1 ^ m0; /* l1 ^= l0 ^ h0 ^ m0; */
}
/* Compute xor-multiply of two binary polynomials (a2, a1, a0) x (b2, b1, b0)
* result is a binary polynomial in 6 mp_digits r[6].
* The caller MUST ensure that r has the right amount of space allocated.
*/
void
s_bmul_3x3(mp_digit *r, const mp_digit a2, const mp_digit a1, const mp_digit a0,
const mp_digit b2, const mp_digit b1, const mp_digit b0)
{
mp_digit zm[4];
s_bmul_1x1(r + 5, r + 4, a2, b2); /* fill top 2 words */
s_bmul_2x2(zm, a1, a2 ^ a0, b1, b2 ^ b0); /* fill middle 4 words */
s_bmul_2x2(r, a1, a0, b1, b0); /* fill bottom 4 words */
zm[3] ^= r[3];
zm[2] ^= r[2];
zm[1] ^= r[1] ^ r[5];
zm[0] ^= r[0] ^ r[4];
r[5] ^= zm[3];
r[4] ^= zm[2];
r[3] ^= zm[1];
r[2] ^= zm[0];
}
/* Compute xor-multiply of two binary polynomials (a3, a2, a1, a0) x (b3, b2, b1, b0)
* result is a binary polynomial in 8 mp_digits r[8].
* The caller MUST ensure that r has the right amount of space allocated.
*/
void
s_bmul_4x4(mp_digit *r, const mp_digit a3, const mp_digit a2, const mp_digit a1,
const mp_digit a0, const mp_digit b3, const mp_digit b2, const mp_digit b1,
const mp_digit b0)
{
mp_digit zm[4];
s_bmul_2x2(r + 4, a3, a2, b3, b2); /* fill top 4 words */
s_bmul_2x2(zm, a3 ^ a1, a2 ^ a0, b3 ^ b1, b2 ^ b0); /* fill middle 4 words */
s_bmul_2x2(r, a1, a0, b1, b0); /* fill bottom 4 words */
zm[3] ^= r[3] ^ r[7];
zm[2] ^= r[2] ^ r[6];
zm[1] ^= r[1] ^ r[5];
zm[0] ^= r[0] ^ r[4];
r[5] ^= zm[3];
r[4] ^= zm[2];
r[3] ^= zm[1];
r[2] ^= zm[0];
}
/* Compute addition of two binary polynomials a and b,
* store result in c; c could be a or b, a and b could be equal;
* c is the bitwise XOR of a and b.
*/
mp_err
mp_badd(const mp_int *a, const mp_int *b, mp_int *c)
{
mp_digit *pa, *pb, *pc;
mp_size ix;
mp_size used_pa, used_pb;
mp_err res = MP_OKAY;
/* Add all digits up to the precision of b. If b had more
* precision than a initially, swap a, b first
*/
if (MP_USED(a) >= MP_USED(b)) {
pa = MP_DIGITS(a);
pb = MP_DIGITS(b);
used_pa = MP_USED(a);
used_pb = MP_USED(b);
} else {
pa = MP_DIGITS(b);
pb = MP_DIGITS(a);
used_pa = MP_USED(b);
used_pb = MP_USED(a);
}
/* Make sure c has enough precision for the output value */
MP_CHECKOK(s_mp_pad(c, used_pa));
/* Do word-by-word xor */
pc = MP_DIGITS(c);
for (ix = 0; ix < used_pb; ix++) {
(*pc++) = (*pa++) ^ (*pb++);
}
/* Finish the rest of digits until we're actually done */
for (; ix < used_pa; ++ix) {
*pc++ = *pa++;
}
MP_USED(c) = used_pa;
MP_SIGN(c) = ZPOS;
s_mp_clamp(c);
CLEANUP:
return res;
}
#define s_mp_div2(a) MP_CHECKOK(mpl_rsh((a), (a), 1));
/* Compute binary polynomial multiply d = a * b */
static void
s_bmul_d(const mp_digit *a, mp_size a_len, mp_digit b, mp_digit *d)
{
mp_digit a_i, a0b0, a1b1, carry = 0;
while (a_len--) {
a_i = *a++;
s_bmul_1x1(&a1b1, &a0b0, a_i, b);
*d++ = a0b0 ^ carry;
carry = a1b1;
}
*d = carry;
}
/* Compute binary polynomial xor multiply accumulate d ^= a * b */
static void
s_bmul_d_add(const mp_digit *a, mp_size a_len, mp_digit b, mp_digit *d)
{
mp_digit a_i, a0b0, a1b1, carry = 0;
while (a_len--) {
a_i = *a++;
s_bmul_1x1(&a1b1, &a0b0, a_i, b);
*d++ ^= a0b0 ^ carry;
carry = a1b1;
}
*d ^= carry;
}
/* Compute binary polynomial xor multiply c = a * b.
* All parameters may be identical.
*/
mp_err
mp_bmul(const mp_int *a, const mp_int *b, mp_int *c)
{
mp_digit *pb, b_i;
mp_int tmp;
mp_size ib, a_used, b_used;
mp_err res = MP_OKAY;
MP_DIGITS(&tmp) = 0;
ARGCHK(a != NULL && b != NULL && c != NULL, MP_BADARG);
if (a == c) {
MP_CHECKOK(mp_init_copy(&tmp, a));
if (a == b)
b = &tmp;
a = &tmp;
} else if (b == c) {
MP_CHECKOK(mp_init_copy(&tmp, b));
b = &tmp;
}
if (MP_USED(a) < MP_USED(b)) {
const mp_int *xch = b; /* switch a and b if b longer */
b = a;
a = xch;
}
MP_USED(c) = 1;
MP_DIGIT(c, 0) = 0;
MP_CHECKOK(s_mp_pad(c, USED(a) + USED(b)));
pb = MP_DIGITS(b);
s_bmul_d(MP_DIGITS(a), MP_USED(a), *pb++, MP_DIGITS(c));
/* Outer loop: Digits of b */
a_used = MP_USED(a);
b_used = MP_USED(b);
MP_USED(c) = a_used + b_used;
for (ib = 1; ib < b_used; ib++) {
b_i = *pb++;
/* Inner product: Digits of a */
if (b_i)
s_bmul_d_add(MP_DIGITS(a), a_used, b_i, MP_DIGITS(c) + ib);
else
MP_DIGIT(c, ib + a_used) = b_i;
}
s_mp_clamp(c);
SIGN(c) = ZPOS;
CLEANUP:
mp_clear(&tmp);
return res;
}
/* Compute modular reduction of a and store result in r.
* r could be a.
* For modular arithmetic, the irreducible polynomial f(t) is represented
* as an array of int[], where f(t) is of the form:
* f(t) = t^p[0] + t^p[1] + ... + t^p[k]
* where m = p[0] > p[1] > ... > p[k] = 0.
*/
mp_err
mp_bmod(const mp_int *a, const unsigned int p[], mp_int *r)
{
int j, k;
int n, dN, d0, d1;
mp_digit zz, *z, tmp;
mp_size used;
mp_err res = MP_OKAY;
/* The algorithm does the reduction in place in r,
* if a != r, copy a into r first so reduction can be done in r
*/
if (a != r) {
MP_CHECKOK(mp_copy(a, r));
}
z = MP_DIGITS(r);
/* start reduction */
/*dN = p[0] / MP_DIGIT_BITS; */
dN = p[0] >> MP_DIGIT_BITS_LOG_2;
used = MP_USED(r);
for (j = used - 1; j > dN;) {
zz = z[j];
if (zz == 0) {
j--;
continue;
}
z[j] = 0;
for (k = 1; p[k] > 0; k++) {
/* reducing component t^p[k] */
n = p[0] - p[k];
/*d0 = n % MP_DIGIT_BITS; */
d0 = n & MP_DIGIT_BITS_MASK;
d1 = MP_DIGIT_BITS - d0;
/*n /= MP_DIGIT_BITS; */
n >>= MP_DIGIT_BITS_LOG_2;
z[j - n] ^= (zz >> d0);
if (d0)
z[j - n - 1] ^= (zz << d1);
}
/* reducing component t^0 */
n = dN;
/*d0 = p[0] % MP_DIGIT_BITS;*/
d0 = p[0] & MP_DIGIT_BITS_MASK;
d1 = MP_DIGIT_BITS - d0;
z[j - n] ^= (zz >> d0);
if (d0)
z[j - n - 1] ^= (zz << d1);
}
/* final round of reduction */
while (j == dN) {
/* d0 = p[0] % MP_DIGIT_BITS; */
d0 = p[0] & MP_DIGIT_BITS_MASK;
zz = z[dN] >> d0;
if (zz == 0)
break;
d1 = MP_DIGIT_BITS - d0;
/* clear up the top d1 bits */
if (d0) {
z[dN] = (z[dN] << d1) >> d1;
} else {
z[dN] = 0;
}
*z ^= zz; /* reduction t^0 component */
for (k = 1; p[k] > 0; k++) {
/* reducing component t^p[k]*/
/* n = p[k] / MP_DIGIT_BITS; */
n = p[k] >> MP_DIGIT_BITS_LOG_2;
/* d0 = p[k] % MP_DIGIT_BITS; */
d0 = p[k] & MP_DIGIT_BITS_MASK;
d1 = MP_DIGIT_BITS - d0;
z[n] ^= (zz << d0);
tmp = zz >> d1;
if (d0 && tmp)
z[n + 1] ^= tmp;
}
}
s_mp_clamp(r);
CLEANUP:
return res;
}
/* Compute the product of two polynomials a and b, reduce modulo p,
* Store the result in r. r could be a or b; a could be b.
*/
mp_err
mp_bmulmod(const mp_int *a, const mp_int *b, const unsigned int p[], mp_int *r)
{
mp_err res;
if (a == b)
return mp_bsqrmod(a, p, r);
if ((res = mp_bmul(a, b, r)) != MP_OKAY)
return res;
return mp_bmod(r, p, r);
}
/* Compute binary polynomial squaring c = a*a mod p .
* Parameter r and a can be identical.
*/
mp_err
mp_bsqrmod(const mp_int *a, const unsigned int p[], mp_int *r)
{
mp_digit *pa, *pr, a_i;
mp_int tmp;
mp_size ia, a_used;
mp_err res;
ARGCHK(a != NULL && r != NULL, MP_BADARG);
MP_DIGITS(&tmp) = 0;
if (a == r) {
MP_CHECKOK(mp_init_copy(&tmp, a));
a = &tmp;
}
MP_USED(r) = 1;
MP_DIGIT(r, 0) = 0;
MP_CHECKOK(s_mp_pad(r, 2 * USED(a)));
pa = MP_DIGITS(a);
pr = MP_DIGITS(r);
a_used = MP_USED(a);
MP_USED(r) = 2 * a_used;
for (ia = 0; ia < a_used; ia++) {
a_i = *pa++;
*pr++ = gf2m_SQR0(a_i);
*pr++ = gf2m_SQR1(a_i);
}
MP_CHECKOK(mp_bmod(r, p, r));
s_mp_clamp(r);
SIGN(r) = ZPOS;
CLEANUP:
mp_clear(&tmp);
return res;
}
/* Compute binary polynomial y/x mod p, y divided by x, reduce modulo p.
* Store the result in r. r could be x or y, and x could equal y.
* Uses algorithm Modular_Division_GF(2^m) from
* Chang-Shantz, S. "From Euclid's GCD to Montgomery Multiplication to
* the Great Divide".
*/
int
mp_bdivmod(const mp_int *y, const mp_int *x, const mp_int *pp,
const unsigned int p[], mp_int *r)
{
mp_int aa, bb, uu;
mp_int *a, *b, *u, *v;
mp_err res = MP_OKAY;
MP_DIGITS(&aa) = 0;
MP_DIGITS(&bb) = 0;
MP_DIGITS(&uu) = 0;
MP_CHECKOK(mp_init_copy(&aa, x));
MP_CHECKOK(mp_init_copy(&uu, y));
MP_CHECKOK(mp_init_copy(&bb, pp));
MP_CHECKOK(s_mp_pad(r, USED(pp)));
MP_USED(r) = 1;
MP_DIGIT(r, 0) = 0;
a = &aa;
b = &bb;
u = &uu;
v = r;
/* reduce x and y mod p */
MP_CHECKOK(mp_bmod(a, p, a));
MP_CHECKOK(mp_bmod(u, p, u));
while (!mp_isodd(a)) {
s_mp_div2(a);
if (mp_isodd(u)) {
MP_CHECKOK(mp_badd(u, pp, u));
}
s_mp_div2(u);
}
do {
if (mp_cmp_mag(b, a) > 0) {
MP_CHECKOK(mp_badd(b, a, b));
MP_CHECKOK(mp_badd(v, u, v));
do {
s_mp_div2(b);
if (mp_isodd(v)) {
MP_CHECKOK(mp_badd(v, pp, v));
}
s_mp_div2(v);
} while (!mp_isodd(b));
} else if ((MP_DIGIT(a, 0) == 1) && (MP_USED(a) == 1))
break;
else {
MP_CHECKOK(mp_badd(a, b, a));
MP_CHECKOK(mp_badd(u, v, u));
do {
s_mp_div2(a);
if (mp_isodd(u)) {
MP_CHECKOK(mp_badd(u, pp, u));
}
s_mp_div2(u);
} while (!mp_isodd(a));
}
} while (1);
MP_CHECKOK(mp_copy(u, r));
CLEANUP:
mp_clear(&aa);
mp_clear(&bb);
mp_clear(&uu);
return res;
}
/* Convert the bit-string representation of a polynomial a into an array
* of integers corresponding to the bits with non-zero coefficient.
* Up to max elements of the array will be filled. Return value is total
* number of coefficients that would be extracted if array was large enough.
*/
int
mp_bpoly2arr(const mp_int *a, unsigned int p[], int max)
{
int i, j, k;
mp_digit top_bit, mask;
top_bit = 1;
top_bit <<= MP_DIGIT_BIT - 1;
for (k = 0; k < max; k++)
p[k] = 0;
k = 0;
for (i = MP_USED(a) - 1; i >= 0; i--) {
mask = top_bit;
for (j = MP_DIGIT_BIT - 1; j >= 0; j--) {
if (MP_DIGITS(a)[i] & mask) {
if (k < max)
p[k] = MP_DIGIT_BIT * i + j;
k++;
}
mask >>= 1;
}
}
return k;
}
/* Convert the coefficient array representation of a polynomial to a
* bit-string. The array must be terminated by 0.
*/
mp_err
mp_barr2poly(const unsigned int p[], mp_int *a)
{
mp_err res = MP_OKAY;
int i;
mp_zero(a);
for (i = 0; p[i] > 0; i++) {
MP_CHECKOK(mpl_set_bit(a, p[i], 1));
}
MP_CHECKOK(mpl_set_bit(a, 0, 1));
CLEANUP:
return res;
}
Ссылка в новой задаче Показать git blame Копировать постоянную ссылку