diff options
Diffstat (limited to 'test/monniaux/BearSSL/src/int')
56 files changed, 5289 insertions, 0 deletions
diff --git a/test/monniaux/BearSSL/src/int/i15_add.c b/test/monniaux/BearSSL/src/int/i15_add.c new file mode 100644 index 00000000..97e29b82 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_add.c @@ -0,0 +1,46 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i15_add(uint16_t *a, const uint16_t *b, uint32_t ctl) +{ + uint32_t cc; + size_t u, m; + + cc = 0; + m = (a[0] + 31) >> 4; + for (u = 1; u < m; u ++) { + uint32_t aw, bw, naw; + + aw = a[u]; + bw = b[u]; + naw = aw + bw + cc; + cc = naw >> 15; + a[u] = MUX(ctl, naw & 0x7FFF, aw); + } + return cc; +} diff --git a/test/monniaux/BearSSL/src/int/i15_bitlen.c b/test/monniaux/BearSSL/src/int/i15_bitlen.c new file mode 100644 index 00000000..ad744671 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_bitlen.c @@ -0,0 +1,44 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i15_bit_length(uint16_t *x, size_t xlen) +{ + uint32_t tw, twk; + + tw = 0; + twk = 0; + while (xlen -- > 0) { + uint32_t w, c; + + c = EQ(tw, 0); + w = x[xlen]; + tw = MUX(c, w, tw); + twk = MUX(c, (uint32_t)xlen, twk); + } + return (twk << 4) + BIT_LENGTH(tw); +} diff --git a/test/monniaux/BearSSL/src/int/i15_decmod.c b/test/monniaux/BearSSL/src/int/i15_decmod.c new file mode 100644 index 00000000..6076c57b --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_decmod.c @@ -0,0 +1,124 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i15_decode_mod(uint16_t *x, const void *src, size_t len, const uint16_t *m) +{ + /* + * Two-pass algorithm: in the first pass, we determine whether the + * value fits; in the second pass, we do the actual write. + * + * During the first pass, 'r' contains the comparison result so + * far: + * 0x00000000 value is equal to the modulus + * 0x00000001 value is greater than the modulus + * 0xFFFFFFFF value is lower than the modulus + * + * Since we iterate starting with the least significant bytes (at + * the end of src[]), each new comparison overrides the previous + * except when the comparison yields 0 (equal). + * + * During the second pass, 'r' is either 0xFFFFFFFF (value fits) + * or 0x00000000 (value does not fit). + * + * We must iterate over all bytes of the source, _and_ possibly + * some extra virtual bytes (with value 0) so as to cover the + * complete modulus as well. We also add 4 such extra bytes beyond + * the modulus length because it then guarantees that no accumulated + * partial word remains to be processed. + */ + const unsigned char *buf; + size_t mlen, tlen; + int pass; + uint32_t r; + + buf = src; + mlen = (m[0] + 15) >> 4; + tlen = (mlen << 1); + if (tlen < len) { + tlen = len; + } + tlen += 4; + r = 0; + for (pass = 0; pass < 2; pass ++) { + size_t u, v; + uint32_t acc; + int acc_len; + + v = 1; + acc = 0; + acc_len = 0; + for (u = 0; u < tlen; u ++) { + uint32_t b; + + if (u < len) { + b = buf[len - 1 - u]; + } else { + b = 0; + } + acc |= (b << acc_len); + acc_len += 8; + if (acc_len >= 15) { + uint32_t xw; + + xw = acc & (uint32_t)0x7FFF; + acc_len -= 15; + acc = b >> (8 - acc_len); + if (v <= mlen) { + if (pass) { + x[v] = r & xw; + } else { + uint32_t cc; + + cc = (uint32_t)CMP(xw, m[v]); + r = MUX(EQ(cc, 0), r, cc); + } + } else { + if (!pass) { + r = MUX(EQ(xw, 0), r, 1); + } + } + v ++; + } + } + + /* + * When we reach this point at the end of the first pass: + * r is either 0, 1 or -1; we want to set r to 0 if it + * is equal to 0 or 1, and leave it to -1 otherwise. + * + * When we reach this point at the end of the second pass: + * r is either 0 or -1; we want to leave that value + * untouched. This is a subcase of the previous. + */ + r >>= 1; + r |= (r << 1); + } + + x[0] = m[0]; + return r & (uint32_t)1; +} diff --git a/test/monniaux/BearSSL/src/int/i15_decode.c b/test/monniaux/BearSSL/src/int/i15_decode.c new file mode 100644 index 00000000..fc2c0be0 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_decode.c @@ -0,0 +1,56 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i15_decode(uint16_t *x, const void *src, size_t len) +{ + const unsigned char *buf; + size_t v; + uint32_t acc; + int acc_len; + + buf = src; + v = 1; + acc = 0; + acc_len = 0; + while (len -- > 0) { + uint32_t b; + + b = buf[len]; + acc |= (b << acc_len); + acc_len += 8; + if (acc_len >= 15) { + x[v ++] = acc & 0x7FFF; + acc_len -= 15; + acc >>= 15; + } + } + if (acc_len != 0) { + x[v ++] = acc; + } + x[0] = br_i15_bit_length(x + 1, v - 1); +} diff --git a/test/monniaux/BearSSL/src/int/i15_decred.c b/test/monniaux/BearSSL/src/int/i15_decred.c new file mode 100644 index 00000000..81e7dd1b --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_decred.c @@ -0,0 +1,100 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i15_decode_reduce(uint16_t *x, + const void *src, size_t len, const uint16_t *m) +{ + uint32_t m_ebitlen, m_rbitlen; + size_t mblen, k; + const unsigned char *buf; + uint32_t acc; + int acc_len; + + /* + * Get the encoded bit length. + */ + m_ebitlen = m[0]; + + /* + * Special case for an invalid (null) modulus. + */ + if (m_ebitlen == 0) { + x[0] = 0; + return; + } + + /* + * Clear the destination. + */ + br_i15_zero(x, m_ebitlen); + + /* + * First decode directly as many bytes as possible. This requires + * computing the actual bit length. + */ + m_rbitlen = m_ebitlen >> 4; + m_rbitlen = (m_ebitlen & 15) + (m_rbitlen << 4) - m_rbitlen; + mblen = (m_rbitlen + 7) >> 3; + k = mblen - 1; + if (k >= len) { + br_i15_decode(x, src, len); + x[0] = m_ebitlen; + return; + } + buf = src; + br_i15_decode(x, buf, k); + x[0] = m_ebitlen; + + /* + * Input remaining bytes, using 15-bit words. + */ + acc = 0; + acc_len = 0; + while (k < len) { + uint32_t v; + + v = buf[k ++]; + acc = (acc << 8) | v; + acc_len += 8; + if (acc_len >= 15) { + br_i15_muladd_small(x, acc >> (acc_len - 15), m); + acc_len -= 15; + acc &= ~((uint32_t)-1 << acc_len); + } + } + + /* + * We may have some bits accumulated. We then perform a shift to + * be able to inject these bits as a full 15-bit word. + */ + if (acc_len != 0) { + acc = (acc | (x[1] << acc_len)) & 0x7FFF; + br_i15_rshift(x, 15 - acc_len); + br_i15_muladd_small(x, acc, m); + } +} diff --git a/test/monniaux/BearSSL/src/int/i15_encode.c b/test/monniaux/BearSSL/src/int/i15_encode.c new file mode 100644 index 00000000..50668f47 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_encode.c @@ -0,0 +1,56 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i15_encode(void *dst, size_t len, const uint16_t *x) +{ + unsigned char *buf; + size_t u, xlen; + uint32_t acc; + int acc_len; + + xlen = (x[0] + 15) >> 4; + if (xlen == 0) { + memset(dst, 0, len); + return; + } + u = 1; + acc = 0; + acc_len = 0; + buf = dst; + while (len -- > 0) { + if (acc_len < 8) { + if (u <= xlen) { + acc += (uint32_t)x[u ++] << acc_len; + } + acc_len += 15; + } + buf[len] = (unsigned char)acc; + acc >>= 8; + acc_len -= 8; + } +} diff --git a/test/monniaux/BearSSL/src/int/i15_fmont.c b/test/monniaux/BearSSL/src/int/i15_fmont.c new file mode 100644 index 00000000..3450b723 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_fmont.c @@ -0,0 +1,59 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i15_from_monty(uint16_t *x, const uint16_t *m, uint16_t m0i) +{ + size_t len, u, v; + + len = (m[0] + 15) >> 4; + for (u = 0; u < len; u ++) { + uint32_t f, cc; + + f = MUL15(x[1], m0i) & 0x7FFF; + cc = 0; + for (v = 0; v < len; v ++) { + uint32_t z; + + z = (uint32_t)x[v + 1] + MUL15(f, m[v + 1]) + cc; + cc = z >> 15; + if (v != 0) { + x[v] = z & 0x7FFF; + } + } + x[len] = cc; + } + + /* + * We may have to do an extra subtraction, but only if the + * value in x[] is indeed greater than or equal to that of m[], + * which is why we must do two calls (first call computes the + * carry, second call performs the subtraction only if the carry + * is 0). + */ + br_i15_sub(x, m, NOT(br_i15_sub(x, m, 0))); +} diff --git a/test/monniaux/BearSSL/src/int/i15_iszero.c b/test/monniaux/BearSSL/src/int/i15_iszero.c new file mode 100644 index 00000000..d4b6f10b --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_iszero.c @@ -0,0 +1,39 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i15_iszero(const uint16_t *x) +{ + uint32_t z; + size_t u; + + z = 0; + for (u = (x[0] + 15) >> 4; u > 0; u --) { + z |= x[u]; + } + return ~(z | -z) >> 31; +} diff --git a/test/monniaux/BearSSL/src/int/i15_moddiv.c b/test/monniaux/BearSSL/src/int/i15_moddiv.c new file mode 100644 index 00000000..45af756d --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_moddiv.c @@ -0,0 +1,465 @@ +/* + * Copyright (c) 2018 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* + * In this file, we handle big integers with a custom format, i.e. + * without the usual one-word header. Value is split into 15-bit words, + * each stored in a 16-bit slot (top bit is zero) in little-endian + * order. The length (in words) is provided explicitly. In some cases, + * the value can be negative (using two's complement representation). In + * some cases, the top word is allowed to have a 16th bit. + */ + +/* + * Negate big integer conditionally. The value consists of 'len' words, + * with 15 bits in each word (the top bit of each word should be 0, + * except possibly for the last word). If 'ctl' is 1, the negation is + * computed; otherwise, if 'ctl' is 0, then the value is unchanged. + */ +static void +cond_negate(uint16_t *a, size_t len, uint32_t ctl) +{ + size_t k; + uint32_t cc, xm; + + cc = ctl; + xm = 0x7FFF & -ctl; + for (k = 0; k < len; k ++) { + uint32_t aw; + + aw = a[k]; + aw = (aw ^ xm) + cc; + a[k] = aw & 0x7FFF; + cc = (aw >> 15) & 1; + } +} + +/* + * Finish modular reduction. Rules on input parameters: + * + * if neg = 1, then -m <= a < 0 + * if neg = 0, then 0 <= a < 2*m + * + * If neg = 0, then the top word of a[] may use 16 bits. + * + * Also, modulus m must be odd. + */ +static void +finish_mod(uint16_t *a, size_t len, const uint16_t *m, uint32_t neg) +{ + size_t k; + uint32_t cc, xm, ym; + + /* + * First pass: compare a (assumed nonnegative) with m. + */ + cc = 0; + for (k = 0; k < len; k ++) { + uint32_t aw, mw; + + aw = a[k]; + mw = m[k]; + cc = (aw - mw - cc) >> 31; + } + + /* + * At this point: + * if neg = 1, then we must add m (regardless of cc) + * if neg = 0 and cc = 0, then we must subtract m + * if neg = 0 and cc = 1, then we must do nothing + */ + xm = 0x7FFF & -neg; + ym = -(neg | (1 - cc)); + cc = neg; + for (k = 0; k < len; k ++) { + uint32_t aw, mw; + + aw = a[k]; + mw = (m[k] ^ xm) & ym; + aw = aw - mw - cc; + a[k] = aw & 0x7FFF; + cc = aw >> 31; + } +} + +/* + * Compute: + * a <- (a*pa+b*pb)/(2^15) + * b <- (a*qa+b*qb)/(2^15) + * The division is assumed to be exact (i.e. the low word is dropped). + * If the final a is negative, then it is negated. Similarly for b. + * Returned value is the combination of two bits: + * bit 0: 1 if a had to be negated, 0 otherwise + * bit 1: 1 if b had to be negated, 0 otherwise + * + * Factors pa, pb, qa and qb must be at most 2^15 in absolute value. + * Source integers a and b must be nonnegative; top word is not allowed + * to contain an extra 16th bit. + */ +static uint32_t +co_reduce(uint16_t *a, uint16_t *b, size_t len, + int32_t pa, int32_t pb, int32_t qa, int32_t qb) +{ + size_t k; + int32_t cca, ccb; + uint32_t nega, negb; + + cca = 0; + ccb = 0; + for (k = 0; k < len; k ++) { + uint32_t wa, wb, za, zb; + uint16_t tta, ttb; + + /* + * Since: + * |pa| <= 2^15 + * |pb| <= 2^15 + * 0 <= wa <= 2^15 - 1 + * 0 <= wb <= 2^15 - 1 + * |cca| <= 2^16 - 1 + * Then: + * |za| <= (2^15-1)*(2^16) + (2^16-1) = 2^31 - 1 + * + * Thus, the new value of cca is such that |cca| <= 2^16 - 1. + * The same applies to ccb. + */ + wa = a[k]; + wb = b[k]; + za = wa * (uint32_t)pa + wb * (uint32_t)pb + (uint32_t)cca; + zb = wa * (uint32_t)qa + wb * (uint32_t)qb + (uint32_t)ccb; + if (k > 0) { + a[k - 1] = za & 0x7FFF; + b[k - 1] = zb & 0x7FFF; + } + tta = za >> 15; + ttb = zb >> 15; + cca = *(int16_t *)&tta; + ccb = *(int16_t *)&ttb; + } + a[len - 1] = (uint16_t)cca; + b[len - 1] = (uint16_t)ccb; + nega = (uint32_t)cca >> 31; + negb = (uint32_t)ccb >> 31; + cond_negate(a, len, nega); + cond_negate(b, len, negb); + return nega | (negb << 1); +} + +/* + * Compute: + * a <- (a*pa+b*pb)/(2^15) mod m + * b <- (a*qa+b*qb)/(2^15) mod m + * + * m0i is equal to -1/m[0] mod 2^15. + * + * Factors pa, pb, qa and qb must be at most 2^15 in absolute value. + * Source integers a and b must be nonnegative; top word is not allowed + * to contain an extra 16th bit. + */ +static void +co_reduce_mod(uint16_t *a, uint16_t *b, size_t len, + int32_t pa, int32_t pb, int32_t qa, int32_t qb, + const uint16_t *m, uint16_t m0i) +{ + size_t k; + int32_t cca, ccb, fa, fb; + + cca = 0; + ccb = 0; + fa = ((a[0] * (uint32_t)pa + b[0] * (uint32_t)pb) * m0i) & 0x7FFF; + fb = ((a[0] * (uint32_t)qa + b[0] * (uint32_t)qb) * m0i) & 0x7FFF; + for (k = 0; k < len; k ++) { + uint32_t wa, wb, za, zb; + uint32_t tta, ttb; + + /* + * In this loop, carries 'cca' and 'ccb' always fit on + * 17 bits (in absolute value). + */ + wa = a[k]; + wb = b[k]; + za = wa * (uint32_t)pa + wb * (uint32_t)pb + + m[k] * (uint32_t)fa + (uint32_t)cca; + zb = wa * (uint32_t)qa + wb * (uint32_t)qb + + m[k] * (uint32_t)fb + (uint32_t)ccb; + if (k > 0) { + a[k - 1] = za & 0x7FFF; + b[k - 1] = zb & 0x7FFF; + } + + /* + * The XOR-and-sub construction below does an arithmetic + * right shift in a portable way (technically, right-shifting + * a negative signed value is implementation-defined in C). + */ +#define M ((uint32_t)1 << 16) + tta = za >> 15; + ttb = zb >> 15; + tta = (tta ^ M) - M; + ttb = (ttb ^ M) - M; + cca = *(int32_t *)&tta; + ccb = *(int32_t *)&ttb; +#undef M + } + a[len - 1] = (uint32_t)cca; + b[len - 1] = (uint32_t)ccb; + + /* + * At this point: + * -m <= a < 2*m + * -m <= b < 2*m + * (this is a case of Montgomery reduction) + * The top word of 'a' and 'b' may have a 16-th bit set. + * We may have to add or subtract the modulus. + */ + finish_mod(a, len, m, (uint32_t)cca >> 31); + finish_mod(b, len, m, (uint32_t)ccb >> 31); +} + +/* see inner.h */ +uint32_t +br_i15_moddiv(uint16_t *x, const uint16_t *y, const uint16_t *m, uint16_t m0i, + uint16_t *t) +{ + /* + * Algorithm is an extended binary GCD. We maintain four values + * a, b, u and v, with the following invariants: + * + * a * x = y * u mod m + * b * x = y * v mod m + * + * Starting values are: + * + * a = y + * b = m + * u = x + * v = 0 + * + * The formal definition of the algorithm is a sequence of steps: + * + * - If a is even, then a <- a/2 and u <- u/2 mod m. + * - Otherwise, if b is even, then b <- b/2 and v <- v/2 mod m. + * - Otherwise, if a > b, then a <- (a-b)/2 and u <- (u-v)/2 mod m. + * - Otherwise, b <- (b-a)/2 and v <- (v-u)/2 mod m. + * + * Algorithm stops when a = b. At that point, they both are equal + * to GCD(y,m); the modular division succeeds if that value is 1. + * The result of the modular division is then u (or v: both are + * equal at that point). + * + * Each step makes either a or b shrink by at least one bit; hence, + * if m has bit length k bits, then 2k-2 steps are sufficient. + * + * + * Though complexity is quadratic in the size of m, the bit-by-bit + * processing is not very efficient. We can speed up processing by + * remarking that the decisions are taken based only on observation + * of the top and low bits of a and b. + * + * In the loop below, at each iteration, we use the two top words + * of a and b, and the low words of a and b, to compute reduction + * parameters pa, pb, qa and qb such that the new values for a + * and b are: + * + * a' = (a*pa + b*pb) / (2^15) + * b' = (a*qa + b*qb) / (2^15) + * + * the division being exact. + * + * Since the choices are based on the top words, they may be slightly + * off, requiring an optional correction: if a' < 0, then we replace + * pa with -pa, and pb with -pb. The total length of a and b is + * thus reduced by at least 14 bits at each iteration. + * + * The stopping conditions are still the same, though: when a + * and b become equal, they must be both odd (since m is odd, + * the GCD cannot be even), therefore the next operation is a + * subtraction, and one of the values becomes 0. At that point, + * nothing else happens, i.e. one value is stuck at 0, and the + * other one is the GCD. + */ + size_t len, k; + uint16_t *a, *b, *u, *v; + uint32_t num, r; + + len = (m[0] + 15) >> 4; + a = t; + b = a + len; + u = x + 1; + v = b + len; + memcpy(a, y + 1, len * sizeof *y); + memcpy(b, m + 1, len * sizeof *m); + memset(v, 0, len * sizeof *v); + + /* + * Loop below ensures that a and b are reduced by some bits each, + * for a total of at least 14 bits. + */ + for (num = ((m[0] - (m[0] >> 4)) << 1) + 14; num >= 14; num -= 14) { + size_t j; + uint32_t c0, c1; + uint32_t a0, a1, b0, b1; + uint32_t a_hi, b_hi, a_lo, b_lo; + int32_t pa, pb, qa, qb; + int i; + + /* + * Extract top words of a and b. If j is the highest + * index >= 1 such that a[j] != 0 or b[j] != 0, then we want + * (a[j] << 15) + a[j - 1], and (b[j] << 15) + b[j - 1]. + * If a and b are down to one word each, then we use a[0] + * and b[0]. + */ + c0 = (uint32_t)-1; + c1 = (uint32_t)-1; + a0 = 0; + a1 = 0; + b0 = 0; + b1 = 0; + j = len; + while (j -- > 0) { + uint32_t aw, bw; + + aw = a[j]; + bw = b[j]; + a0 ^= (a0 ^ aw) & c0; + a1 ^= (a1 ^ aw) & c1; + b0 ^= (b0 ^ bw) & c0; + b1 ^= (b1 ^ bw) & c1; + c1 = c0; + c0 &= (((aw | bw) + 0xFFFF) >> 16) - (uint32_t)1; + } + + /* + * If c1 = 0, then we grabbed two words for a and b. + * If c1 != 0 but c0 = 0, then we grabbed one word. It + * is not possible that c1 != 0 and c0 != 0, because that + * would mean that both integers are zero. + */ + a1 |= a0 & c1; + a0 &= ~c1; + b1 |= b0 & c1; + b0 &= ~c1; + a_hi = (a0 << 15) + a1; + b_hi = (b0 << 15) + b1; + a_lo = a[0]; + b_lo = b[0]; + + /* + * Compute reduction factors: + * + * a' = a*pa + b*pb + * b' = a*qa + b*qb + * + * such that a' and b' are both multiple of 2^15, but are + * only marginally larger than a and b. + */ + pa = 1; + pb = 0; + qa = 0; + qb = 1; + for (i = 0; i < 15; i ++) { + /* + * At each iteration: + * + * a <- (a-b)/2 if: a is odd, b is odd, a_hi > b_hi + * b <- (b-a)/2 if: a is odd, b is odd, a_hi <= b_hi + * a <- a/2 if: a is even + * b <- b/2 if: a is odd, b is even + * + * We multiply a_lo and b_lo by 2 at each + * iteration, thus a division by 2 really is a + * non-multiplication by 2. + */ + uint32_t r, oa, ob, cAB, cBA, cA; + + /* + * cAB = 1 if b must be subtracted from a + * cBA = 1 if a must be subtracted from b + * cA = 1 if a is divided by 2, 0 otherwise + * + * Rules: + * + * cAB and cBA cannot be both 1. + * if a is not divided by 2, b is. + */ + r = GT(a_hi, b_hi); + oa = (a_lo >> i) & 1; + ob = (b_lo >> i) & 1; + cAB = oa & ob & r; + cBA = oa & ob & NOT(r); + cA = cAB | NOT(oa); + + /* + * Conditional subtractions. + */ + a_lo -= b_lo & -cAB; + a_hi -= b_hi & -cAB; + pa -= qa & -(int32_t)cAB; + pb -= qb & -(int32_t)cAB; + b_lo -= a_lo & -cBA; + b_hi -= a_hi & -cBA; + qa -= pa & -(int32_t)cBA; + qb -= pb & -(int32_t)cBA; + + /* + * Shifting. + */ + a_lo += a_lo & (cA - 1); + pa += pa & ((int32_t)cA - 1); + pb += pb & ((int32_t)cA - 1); + a_hi ^= (a_hi ^ (a_hi >> 1)) & -cA; + b_lo += b_lo & -cA; + qa += qa & -(int32_t)cA; + qb += qb & -(int32_t)cA; + b_hi ^= (b_hi ^ (b_hi >> 1)) & (cA - 1); + } + + /* + * Replace a and b with new values a' and b'. + */ + r = co_reduce(a, b, len, pa, pb, qa, qb); + pa -= pa * ((r & 1) << 1); + pb -= pb * ((r & 1) << 1); + qa -= qa * (r & 2); + qb -= qb * (r & 2); + co_reduce_mod(u, v, len, pa, pb, qa, qb, m + 1, m0i); + } + + /* + * Now one of the arrays should be 0, and the other contains + * the GCD. If a is 0, then u is 0 as well, and v contains + * the division result. + * Result is correct if and only if GCD is 1. + */ + r = (a[0] | b[0]) ^ 1; + u[0] |= v[0]; + for (k = 1; k < len; k ++) { + r |= a[k] | b[k]; + u[k] |= v[k]; + } + return EQ0(r); +} diff --git a/test/monniaux/BearSSL/src/int/i15_modpow.c b/test/monniaux/BearSSL/src/int/i15_modpow.c new file mode 100644 index 00000000..9bf304e5 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_modpow.c @@ -0,0 +1,50 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i15_modpow(uint16_t *x, + const unsigned char *e, size_t elen, + const uint16_t *m, uint16_t m0i, uint16_t *t1, uint16_t *t2) +{ + size_t mlen; + unsigned k; + + mlen = ((m[0] + 31) >> 4) * sizeof m[0]; + memcpy(t1, x, mlen); + br_i15_to_monty(t1, m); + br_i15_zero(x, m[0]); + x[1] = 1; + for (k = 0; k < ((unsigned)elen << 3); k ++) { + uint32_t ctl; + + ctl = (e[elen - 1 - (k >> 3)] >> (k & 7)) & 1; + br_i15_montymul(t2, x, t1, m, m0i); + CCOPY(ctl, x, t2, mlen); + br_i15_montymul(t2, t1, t1, m, m0i); + memcpy(t1, t2, mlen); + } +} diff --git a/test/monniaux/BearSSL/src/int/i15_modpow2.c b/test/monniaux/BearSSL/src/int/i15_modpow2.c new file mode 100644 index 00000000..4b321186 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_modpow2.c @@ -0,0 +1,160 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i15_modpow_opt(uint16_t *x, + const unsigned char *e, size_t elen, + const uint16_t *m, uint16_t m0i, uint16_t *tmp, size_t twlen) +{ + size_t mlen, mwlen; + uint16_t *t1, *t2, *base; + size_t u, v; + uint32_t acc; + int acc_len, win_len; + + /* + * Get modulus size. + */ + mwlen = (m[0] + 31) >> 4; + mlen = mwlen * sizeof m[0]; + mwlen += (mwlen & 1); + t1 = tmp; + t2 = tmp + mwlen; + + /* + * Compute possible window size, with a maximum of 5 bits. + * When the window has size 1 bit, we use a specific code + * that requires only two temporaries. Otherwise, for a + * window of k bits, we need 2^k+1 temporaries. + */ + if (twlen < (mwlen << 1)) { + return 0; + } + for (win_len = 5; win_len > 1; win_len --) { + if ((((uint32_t)1 << win_len) + 1) * mwlen <= twlen) { + break; + } + } + + /* + * Everything is done in Montgomery representation. + */ + br_i15_to_monty(x, m); + + /* + * Compute window contents. If the window has size one bit only, + * then t2 is set to x; otherwise, t2[0] is left untouched, and + * t2[k] is set to x^k (for k >= 1). + */ + if (win_len == 1) { + memcpy(t2, x, mlen); + } else { + memcpy(t2 + mwlen, x, mlen); + base = t2 + mwlen; + for (u = 2; u < ((unsigned)1 << win_len); u ++) { + br_i15_montymul(base + mwlen, base, x, m, m0i); + base += mwlen; + } + } + + /* + * We need to set x to 1, in Montgomery representation. This can + * be done efficiently by setting the high word to 1, then doing + * one word-sized shift. + */ + br_i15_zero(x, m[0]); + x[(m[0] + 15) >> 4] = 1; + br_i15_muladd_small(x, 0, m); + + /* + * We process bits from most to least significant. At each + * loop iteration, we have acc_len bits in acc. + */ + acc = 0; + acc_len = 0; + while (acc_len > 0 || elen > 0) { + int i, k; + uint32_t bits; + + /* + * Get the next bits. + */ + k = win_len; + if (acc_len < win_len) { + if (elen > 0) { + acc = (acc << 8) | *e ++; + elen --; + acc_len += 8; + } else { + k = acc_len; + } + } + bits = (acc >> (acc_len - k)) & (((uint32_t)1 << k) - 1); + acc_len -= k; + + /* + * We could get exactly k bits. Compute k squarings. + */ + for (i = 0; i < k; i ++) { + br_i15_montymul(t1, x, x, m, m0i); + memcpy(x, t1, mlen); + } + + /* + * Window lookup: we want to set t2 to the window + * lookup value, assuming the bits are non-zero. If + * the window length is 1 bit only, then t2 is + * already set; otherwise, we do a constant-time lookup. + */ + if (win_len > 1) { + br_i15_zero(t2, m[0]); + base = t2 + mwlen; + for (u = 1; u < ((uint32_t)1 << k); u ++) { + uint32_t mask; + + mask = -EQ(u, bits); + for (v = 1; v < mwlen; v ++) { + t2[v] |= mask & base[v]; + } + base += mwlen; + } + } + + /* + * Multiply with the looked-up value. We keep the + * product only if the exponent bits are not all-zero. + */ + br_i15_montymul(t1, x, t2, m, m0i); + CCOPY(NEQ(bits, 0), x, t1, mlen); + } + + /* + * Convert back from Montgomery representation, and exit. + */ + br_i15_from_monty(x, m, m0i); + return 1; +} diff --git a/test/monniaux/BearSSL/src/int/i15_montmul.c b/test/monniaux/BearSSL/src/int/i15_montmul.c new file mode 100644 index 00000000..e98bc32c --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_montmul.c @@ -0,0 +1,184 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i15_montymul(uint16_t *d, const uint16_t *x, const uint16_t *y, + const uint16_t *m, uint16_t m0i) +{ + size_t len, len4, u, v; + uint32_t dh; + + len = (m[0] + 15) >> 4; + len4 = len & ~(size_t)3; + br_i15_zero(d, m[0]); + dh = 0; + for (u = 0; u < len; u ++) { + uint32_t f, xu, r, zh; + + xu = x[u + 1]; + f = MUL15((d[1] + MUL15(x[u + 1], y[1])) & 0x7FFF, m0i) + & 0x7FFF; +#if BR_ARMEL_CORTEXM_GCC + if (len4 != 0) { + uint16_t *limit; + + limit = d + len4; + asm volatile ( +"\n\ + @ carry: r=r2 \n\ + @ multipliers: xu=r3 f=r4 \n\ + @ base registers: d+v=r5 y+v=r6 m+v=r7 \n\ + @ r8 contains 0x7FFF \n\ + @ r9 contains d+len4 \n\ + ldr r0, %[limit] \n\ + ldr r3, %[xu] \n\ + mov r9, r0 \n\ + ldr r4, %[f] \n\ + eor r2, r2 \n\ + ldr r5, %[d] \n\ + sub r1, r2, #1 \n\ + ldr r6, %[y] \n\ + lsr r1, r1, #17 \n\ + ldr r7, %[m] \n\ + mov r8, r1 \n\ +loop%=: \n\ + ldrh r0, [r6, #2] \n\ + ldrh r1, [r7, #2] \n\ + mul r0, r3 \n\ + mul r1, r4 \n\ + add r2, r0, r2 \n\ + ldrh r0, [r5, #2] \n\ + add r2, r1, r2 \n\ + mov r1, r8 \n\ + add r2, r0, r2 \n\ + and r1, r2 \n\ + lsr r2, r2, #15 \n\ + strh r1, [r5, #0] \n\ + \n\ + ldrh r0, [r6, #4] \n\ + ldrh r1, [r7, #4] \n\ + mul r0, r3 \n\ + mul r1, r4 \n\ + add r2, r0, r2 \n\ + ldrh r0, [r5, #4] \n\ + add r2, r1, r2 \n\ + mov r1, r8 \n\ + add r2, r0, r2 \n\ + and r1, r2 \n\ + lsr r2, r2, #15 \n\ + strh r1, [r5, #2] \n\ + \n\ + ldrh r0, [r6, #6] \n\ + ldrh r1, [r7, #6] \n\ + mul r0, r3 \n\ + mul r1, r4 \n\ + add r2, r0, r2 \n\ + ldrh r0, [r5, #6] \n\ + add r2, r1, r2 \n\ + mov r1, r8 \n\ + add r2, r0, r2 \n\ + and r1, r2 \n\ + lsr r2, r2, #15 \n\ + strh r1, [r5, #4] \n\ + \n\ + ldrh r0, [r6, #8] \n\ + ldrh r1, [r7, #8] \n\ + mul r0, r3 \n\ + mul r1, r4 \n\ + add r2, r0, r2 \n\ + ldrh r0, [r5, #8] \n\ + add r2, r1, r2 \n\ + mov r1, r8 \n\ + add r2, r0, r2 \n\ + and r1, r2 \n\ + lsr r2, r2, #15 \n\ + strh r1, [r5, #6] \n\ + \n\ + add r5, r5, #8 \n\ + add r6, r6, #8 \n\ + add r7, r7, #8 \n\ + cmp r5, r9 \n\ + bne loop%= \n\ + \n\ + str r2, %[carry] \n\ +" +: [carry] "=m" (r) +: [xu] "m" (xu), [f] "m" (f), [d] "m" (d), [y] "m" (y), + [m] "m" (m), [limit] "m" (limit) +: "r0", "r1", "r2", "r3", "r4", "r5", "r6", "r7", "r8", "r9" ); + } else { + r = 0; + } + v = len4; +#else + r = 0; + for (v = 0; v < len4; v += 4) { + uint32_t z; + + z = d[v + 1] + MUL15(xu, y[v + 1]) + + MUL15(f, m[v + 1]) + r; + r = z >> 15; + d[v + 0] = z & 0x7FFF; + z = d[v + 2] + MUL15(xu, y[v + 2]) + + MUL15(f, m[v + 2]) + r; + r = z >> 15; + d[v + 1] = z & 0x7FFF; + z = d[v + 3] + MUL15(xu, y[v + 3]) + + MUL15(f, m[v + 3]) + r; + r = z >> 15; + d[v + 2] = z & 0x7FFF; + z = d[v + 4] + MUL15(xu, y[v + 4]) + + MUL15(f, m[v + 4]) + r; + r = z >> 15; + d[v + 3] = z & 0x7FFF; + } +#endif + for (; v < len; v ++) { + uint32_t z; + + z = d[v + 1] + MUL15(xu, y[v + 1]) + + MUL15(f, m[v + 1]) + r; + r = z >> 15; + d[v + 0] = z & 0x7FFF; + } + + zh = dh + r; + d[len] = zh & 0x7FFF; + dh = zh >> 15; + } + + /* + * Restore the bit length (it was overwritten in the loop above). + */ + d[0] = m[0]; + + /* + * d[] may be greater than m[], but it is still lower than twice + * the modulus. + */ + br_i15_sub(d, m, NEQ(dh, 0) | NOT(br_i15_sub(d, m, 0))); +} diff --git a/test/monniaux/BearSSL/src/int/i15_mulacc.c b/test/monniaux/BearSSL/src/int/i15_mulacc.c new file mode 100644 index 00000000..7a073ac6 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_mulacc.c @@ -0,0 +1,61 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i15_mulacc(uint16_t *d, const uint16_t *a, const uint16_t *b) +{ + size_t alen, blen, u; + unsigned dl, dh; + + alen = (a[0] + 15) >> 4; + blen = (b[0] + 15) >> 4; + + /* + * Announced bit length of d[] will be the sum of the announced + * bit lengths of a[] and b[]; but the lengths are encoded. + */ + dl = (a[0] & 15) + (b[0] & 15); + dh = (a[0] >> 4) + (b[0] >> 4); + d[0] = (dh << 4) + dl + (~(uint32_t)(dl - 15) >> 31); + + for (u = 0; u < blen; u ++) { + uint32_t f; + size_t v; + uint32_t cc; + + f = b[1 + u]; + cc = 0; + for (v = 0; v < alen; v ++) { + uint32_t z; + + z = (uint32_t)d[1 + u + v] + MUL15(f, a[1 + v]) + cc; + cc = z >> 15; + d[1 + u + v] = z & 0x7FFF; + } + d[1 + u + alen] = cc; + } +} diff --git a/test/monniaux/BearSSL/src/int/i15_muladd.c b/test/monniaux/BearSSL/src/int/i15_muladd.c new file mode 100644 index 00000000..c4b72168 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_muladd.c @@ -0,0 +1,173 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* + * Constant-time division. The divisor must not be larger than 16 bits, + * and the quotient must fit on 17 bits. + */ +static uint32_t +divrem16(uint32_t x, uint32_t d, uint32_t *r) +{ + int i; + uint32_t q; + + q = 0; + d <<= 16; + for (i = 16; i >= 0; i --) { + uint32_t ctl; + + ctl = LE(d, x); + q |= ctl << i; + x -= (-ctl) & d; + d >>= 1; + } + if (r != NULL) { + *r = x; + } + return q; +} + +/* see inner.h */ +void +br_i15_muladd_small(uint16_t *x, uint16_t z, const uint16_t *m) +{ + /* + * Constant-time: we accept to leak the exact bit length of the + * modulus m. + */ + unsigned m_bitlen, mblr; + size_t u, mlen; + uint32_t hi, a0, a, b, q; + uint32_t cc, tb, over, under; + + /* + * Simple case: the modulus fits on one word. + */ + m_bitlen = m[0]; + if (m_bitlen == 0) { + return; + } + if (m_bitlen <= 15) { + uint32_t rem; + + divrem16(((uint32_t)x[1] << 15) | z, m[1], &rem); + x[1] = rem; + return; + } + mlen = (m_bitlen + 15) >> 4; + mblr = m_bitlen & 15; + + /* + * Principle: we estimate the quotient (x*2^15+z)/m by + * doing a 30/15 division with the high words. + * + * Let: + * w = 2^15 + * a = (w*a0 + a1) * w^N + a2 + * b = b0 * w^N + b2 + * such that: + * 0 <= a0 < w + * 0 <= a1 < w + * 0 <= a2 < w^N + * w/2 <= b0 < w + * 0 <= b2 < w^N + * a < w*b + * I.e. the two top words of a are a0:a1, the top word of b is + * b0, we ensured that b0 is "full" (high bit set), and a is + * such that the quotient q = a/b fits on one word (0 <= q < w). + * + * If a = b*q + r (with 0 <= r < q), then we can estimate q by + * using a division on the top words: + * a0*w + a1 = b0*u + v (with 0 <= v < b0) + * Then the following holds: + * 0 <= u <= w + * u-2 <= q <= u + */ + hi = x[mlen]; + if (mblr == 0) { + a0 = x[mlen]; + memmove(x + 2, x + 1, (mlen - 1) * sizeof *x); + x[1] = z; + a = (a0 << 15) + x[mlen]; + b = m[mlen]; + } else { + a0 = (x[mlen] << (15 - mblr)) | (x[mlen - 1] >> mblr); + memmove(x + 2, x + 1, (mlen - 1) * sizeof *x); + x[1] = z; + a = (a0 << 15) | (((x[mlen] << (15 - mblr)) + | (x[mlen - 1] >> mblr)) & 0x7FFF); + b = (m[mlen] << (15 - mblr)) | (m[mlen - 1] >> mblr); + } + q = divrem16(a, b, NULL); + + /* + * We computed an estimate for q, but the real one may be q, + * q-1 or q-2; moreover, the division may have returned a value + * 8000 or even 8001 if the two high words were identical, and + * we want to avoid values beyond 7FFF. We thus adjust q so + * that the "true" multiplier will be q+1, q or q-1, and q is + * in the 0000..7FFF range. + */ + q = MUX(EQ(b, a0), 0x7FFF, q - 1 + ((q - 1) >> 31)); + + /* + * We subtract q*m from x (x has an extra high word of value 'hi'). + * Since q may be off by 1 (in either direction), we may have to + * add or subtract m afterwards. + * + * The 'tb' flag will be true (1) at the end of the loop if the + * result is greater than or equal to the modulus (not counting + * 'hi' or the carry). + */ + cc = 0; + tb = 1; + for (u = 1; u <= mlen; u ++) { + uint32_t mw, zl, xw, nxw; + + mw = m[u]; + zl = MUL15(mw, q) + cc; + cc = zl >> 15; + zl &= 0x7FFF; + xw = x[u]; + nxw = xw - zl; + cc += nxw >> 31; + nxw &= 0x7FFF; + x[u] = nxw; + tb = MUX(EQ(nxw, mw), tb, GT(nxw, mw)); + } + + /* + * If we underestimated q, then either cc < hi (one extra bit + * beyond the top array word), or cc == hi and tb is true (no + * extra bit, but the result is not lower than the modulus). + * + * If we overestimated q, then cc > hi. + */ + over = GT(cc, hi); + under = ~over & (tb | LT(cc, hi)); + br_i15_add(x, m, over); + br_i15_sub(x, m, under); +} diff --git a/test/monniaux/BearSSL/src/int/i15_ninv15.c b/test/monniaux/BearSSL/src/int/i15_ninv15.c new file mode 100644 index 00000000..de3a3ba8 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_ninv15.c @@ -0,0 +1,38 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint16_t +br_i15_ninv15(uint16_t x) +{ + uint32_t y; + + y = 2 - x; + y = MUL15(y, 2 - MUL15(x, y)); + y = MUL15(y, 2 - MUL15(x, y)); + y = MUL15(y, 2 - MUL15(x, y)); + return MUX(x & 1, -y, 0) & 0x7FFF; +} diff --git a/test/monniaux/BearSSL/src/int/i15_reduce.c b/test/monniaux/BearSSL/src/int/i15_reduce.c new file mode 100644 index 00000000..0931b104 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_reduce.c @@ -0,0 +1,66 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i15_reduce(uint16_t *x, const uint16_t *a, const uint16_t *m) +{ + uint32_t m_bitlen, a_bitlen; + size_t mlen, alen, u; + + m_bitlen = m[0]; + mlen = (m_bitlen + 15) >> 4; + + x[0] = m_bitlen; + if (m_bitlen == 0) { + return; + } + + /* + * If the source is shorter, then simply copy all words from a[] + * and zero out the upper words. + */ + a_bitlen = a[0]; + alen = (a_bitlen + 15) >> 4; + if (a_bitlen < m_bitlen) { + memcpy(x + 1, a + 1, alen * sizeof *a); + for (u = alen; u < mlen; u ++) { + x[u + 1] = 0; + } + return; + } + + /* + * The source length is at least equal to that of the modulus. + * We must thus copy N-1 words, and input the remaining words + * one by one. + */ + memcpy(x + 1, a + 2 + (alen - mlen), (mlen - 1) * sizeof *a); + x[mlen] = 0; + for (u = 1 + alen - mlen; u > 0; u --) { + br_i15_muladd_small(x, a[u], m); + } +} diff --git a/test/monniaux/BearSSL/src/int/i15_rshift.c b/test/monniaux/BearSSL/src/int/i15_rshift.c new file mode 100644 index 00000000..f9991ab6 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_rshift.c @@ -0,0 +1,47 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i15_rshift(uint16_t *x, int count) +{ + size_t u, len; + unsigned r; + + len = (x[0] + 15) >> 4; + if (len == 0) { + return; + } + r = x[1] >> count; + for (u = 2; u <= len; u ++) { + unsigned w; + + w = x[u]; + x[u - 1] = ((w << (15 - count)) | r) & 0x7FFF; + r = w >> count; + } + x[len] = r; +} diff --git a/test/monniaux/BearSSL/src/int/i15_sub.c b/test/monniaux/BearSSL/src/int/i15_sub.c new file mode 100644 index 00000000..1983c4dc --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_sub.c @@ -0,0 +1,46 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i15_sub(uint16_t *a, const uint16_t *b, uint32_t ctl) +{ + uint32_t cc; + size_t u, m; + + cc = 0; + m = (a[0] + 31) >> 4; + for (u = 1; u < m; u ++) { + uint32_t aw, bw, naw; + + aw = a[u]; + bw = b[u]; + naw = aw - bw - cc; + cc = naw >> 31; + a[u] = MUX(ctl, naw & 0x7FFF, aw); + } + return cc; +} diff --git a/test/monniaux/BearSSL/src/int/i15_tmont.c b/test/monniaux/BearSSL/src/int/i15_tmont.c new file mode 100644 index 00000000..d5c4b8b7 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i15_tmont.c @@ -0,0 +1,36 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i15_to_monty(uint16_t *x, const uint16_t *m) +{ + unsigned k; + + for (k = (m[0] + 15) >> 4; k > 0; k --) { + br_i15_muladd_small(x, 0, m); + } +} diff --git a/test/monniaux/BearSSL/src/int/i31_add.c b/test/monniaux/BearSSL/src/int/i31_add.c new file mode 100644 index 00000000..2ca47c6b --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_add.c @@ -0,0 +1,46 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i31_add(uint32_t *a, const uint32_t *b, uint32_t ctl) +{ + uint32_t cc; + size_t u, m; + + cc = 0; + m = (a[0] + 63) >> 5; + for (u = 1; u < m; u ++) { + uint32_t aw, bw, naw; + + aw = a[u]; + bw = b[u]; + naw = aw + bw + cc; + cc = naw >> 31; + a[u] = MUX(ctl, naw & (uint32_t)0x7FFFFFFF, aw); + } + return cc; +} diff --git a/test/monniaux/BearSSL/src/int/i31_bitlen.c b/test/monniaux/BearSSL/src/int/i31_bitlen.c new file mode 100644 index 00000000..3e127c2c --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_bitlen.c @@ -0,0 +1,44 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i31_bit_length(uint32_t *x, size_t xlen) +{ + uint32_t tw, twk; + + tw = 0; + twk = 0; + while (xlen -- > 0) { + uint32_t w, c; + + c = EQ(tw, 0); + w = x[xlen]; + tw = MUX(c, w, tw); + twk = MUX(c, (uint32_t)xlen, twk); + } + return (twk << 5) + BIT_LENGTH(tw); +} diff --git a/test/monniaux/BearSSL/src/int/i31_decmod.c b/test/monniaux/BearSSL/src/int/i31_decmod.c new file mode 100644 index 00000000..3cd7bfe9 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_decmod.c @@ -0,0 +1,124 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i31_decode_mod(uint32_t *x, const void *src, size_t len, const uint32_t *m) +{ + /* + * Two-pass algorithm: in the first pass, we determine whether the + * value fits; in the second pass, we do the actual write. + * + * During the first pass, 'r' contains the comparison result so + * far: + * 0x00000000 value is equal to the modulus + * 0x00000001 value is greater than the modulus + * 0xFFFFFFFF value is lower than the modulus + * + * Since we iterate starting with the least significant bytes (at + * the end of src[]), each new comparison overrides the previous + * except when the comparison yields 0 (equal). + * + * During the second pass, 'r' is either 0xFFFFFFFF (value fits) + * or 0x00000000 (value does not fit). + * + * We must iterate over all bytes of the source, _and_ possibly + * some extra virtual bytes (with value 0) so as to cover the + * complete modulus as well. We also add 4 such extra bytes beyond + * the modulus length because it then guarantees that no accumulated + * partial word remains to be processed. + */ + const unsigned char *buf; + size_t mlen, tlen; + int pass; + uint32_t r; + + buf = src; + mlen = (m[0] + 31) >> 5; + tlen = (mlen << 2); + if (tlen < len) { + tlen = len; + } + tlen += 4; + r = 0; + for (pass = 0; pass < 2; pass ++) { + size_t u, v; + uint32_t acc; + int acc_len; + + v = 1; + acc = 0; + acc_len = 0; + for (u = 0; u < tlen; u ++) { + uint32_t b; + + if (u < len) { + b = buf[len - 1 - u]; + } else { + b = 0; + } + acc |= (b << acc_len); + acc_len += 8; + if (acc_len >= 31) { + uint32_t xw; + + xw = acc & (uint32_t)0x7FFFFFFF; + acc_len -= 31; + acc = b >> (8 - acc_len); + if (v <= mlen) { + if (pass) { + x[v] = r & xw; + } else { + uint32_t cc; + + cc = (uint32_t)CMP(xw, m[v]); + r = MUX(EQ(cc, 0), r, cc); + } + } else { + if (!pass) { + r = MUX(EQ(xw, 0), r, 1); + } + } + v ++; + } + } + + /* + * When we reach this point at the end of the first pass: + * r is either 0, 1 or -1; we want to set r to 0 if it + * is equal to 0 or 1, and leave it to -1 otherwise. + * + * When we reach this point at the end of the second pass: + * r is either 0 or -1; we want to leave that value + * untouched. This is a subcase of the previous. + */ + r >>= 1; + r |= (r << 1); + } + + x[0] = m[0]; + return r & (uint32_t)1; +} diff --git a/test/monniaux/BearSSL/src/int/i31_decode.c b/test/monniaux/BearSSL/src/int/i31_decode.c new file mode 100644 index 00000000..8ec6d908 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_decode.c @@ -0,0 +1,57 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i31_decode(uint32_t *x, const void *src, size_t len) +{ + const unsigned char *buf; + size_t u, v; + uint32_t acc; + int acc_len; + + buf = src; + u = len; + v = 1; + acc = 0; + acc_len = 0; + while (u -- > 0) { + uint32_t b; + + b = buf[u]; + acc |= (b << acc_len); + acc_len += 8; + if (acc_len >= 31) { + x[v ++] = acc & (uint32_t)0x7FFFFFFF; + acc_len -= 31; + acc = b >> (8 - acc_len); + } + } + if (acc_len != 0) { + x[v ++] = acc; + } + x[0] = br_i31_bit_length(x + 1, v - 1); +} diff --git a/test/monniaux/BearSSL/src/int/i31_decred.c b/test/monniaux/BearSSL/src/int/i31_decred.c new file mode 100644 index 00000000..43db6624 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_decred.c @@ -0,0 +1,103 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i31_decode_reduce(uint32_t *x, + const void *src, size_t len, const uint32_t *m) +{ + uint32_t m_ebitlen, m_rbitlen; + size_t mblen, k; + const unsigned char *buf; + uint32_t acc; + int acc_len; + + /* + * Get the encoded bit length. + */ + m_ebitlen = m[0]; + + /* + * Special case for an invalid (null) modulus. + */ + if (m_ebitlen == 0) { + x[0] = 0; + return; + } + + /* + * Clear the destination. + */ + br_i31_zero(x, m_ebitlen); + + /* + * First decode directly as many bytes as possible. This requires + * computing the actual bit length. + */ + m_rbitlen = m_ebitlen >> 5; + m_rbitlen = (m_ebitlen & 31) + (m_rbitlen << 5) - m_rbitlen; + mblen = (m_rbitlen + 7) >> 3; + k = mblen - 1; + if (k >= len) { + br_i31_decode(x, src, len); + x[0] = m_ebitlen; + return; + } + buf = src; + br_i31_decode(x, buf, k); + x[0] = m_ebitlen; + + /* + * Input remaining bytes, using 31-bit words. + */ + acc = 0; + acc_len = 0; + while (k < len) { + uint32_t v; + + v = buf[k ++]; + if (acc_len >= 23) { + acc_len -= 23; + acc <<= (8 - acc_len); + acc |= v >> acc_len; + br_i31_muladd_small(x, acc, m); + acc = v & (0xFF >> (8 - acc_len)); + } else { + acc = (acc << 8) | v; + acc_len += 8; + } + } + + /* + * We may have some bits accumulated. We then perform a shift to + * be able to inject these bits as a full 31-bit word. + */ + if (acc_len != 0) { + acc = (acc | (x[1] << acc_len)) & 0x7FFFFFFF; + br_i31_rshift(x, 31 - acc_len); + br_i31_muladd_small(x, acc, m); + } +} diff --git a/test/monniaux/BearSSL/src/int/i31_encode.c b/test/monniaux/BearSSL/src/int/i31_encode.c new file mode 100644 index 00000000..b6b40c45 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_encode.c @@ -0,0 +1,79 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i31_encode(void *dst, size_t len, const uint32_t *x) +{ + unsigned char *buf; + size_t k, xlen; + uint32_t acc; + int acc_len; + + xlen = (x[0] + 31) >> 5; + if (xlen == 0) { + memset(dst, 0, len); + return; + } + buf = (unsigned char *)dst + len; + k = 1; + acc = 0; + acc_len = 0; + while (len != 0) { + uint32_t w; + + w = (k <= xlen) ? x[k] : 0; + k ++; + if (acc_len == 0) { + acc = w; + acc_len = 31; + } else { + uint32_t z; + + z = acc | (w << acc_len); + acc_len --; + acc = w >> (31 - acc_len); + if (len >= 4) { + buf -= 4; + len -= 4; + br_enc32be(buf, z); + } else { + switch (len) { + case 3: + buf[-3] = (unsigned char)(z >> 16); + /* fall through */ + case 2: + buf[-2] = (unsigned char)(z >> 8); + /* fall through */ + case 1: + buf[-1] = (unsigned char)z; + break; + } + return; + } + } + } +} diff --git a/test/monniaux/BearSSL/src/int/i31_fmont.c b/test/monniaux/BearSSL/src/int/i31_fmont.c new file mode 100644 index 00000000..c24b4176 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_fmont.c @@ -0,0 +1,60 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i31_from_monty(uint32_t *x, const uint32_t *m, uint32_t m0i) +{ + size_t len, u, v; + + len = (m[0] + 31) >> 5; + for (u = 0; u < len; u ++) { + uint32_t f; + uint64_t cc; + + f = MUL31_lo(x[1], m0i); + cc = 0; + for (v = 0; v < len; v ++) { + uint64_t z; + + z = (uint64_t)x[v + 1] + MUL31(f, m[v + 1]) + cc; + cc = z >> 31; + if (v != 0) { + x[v] = (uint32_t)z & 0x7FFFFFFF; + } + } + x[len] = (uint32_t)cc; + } + + /* + * We may have to do an extra subtraction, but only if the + * value in x[] is indeed greater than or equal to that of m[], + * which is why we must do two calls (first call computes the + * carry, second call performs the subtraction only if the carry + * is 0). + */ + br_i31_sub(x, m, NOT(br_i31_sub(x, m, 0))); +} diff --git a/test/monniaux/BearSSL/src/int/i31_iszero.c b/test/monniaux/BearSSL/src/int/i31_iszero.c new file mode 100644 index 00000000..8a7ea44f --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_iszero.c @@ -0,0 +1,39 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i31_iszero(const uint32_t *x) +{ + uint32_t z; + size_t u; + + z = 0; + for (u = (x[0] + 31) >> 5; u > 0; u --) { + z |= x[u]; + } + return ~(z | -z) >> 31; +} diff --git a/test/monniaux/BearSSL/src/int/i31_moddiv.c b/test/monniaux/BearSSL/src/int/i31_moddiv.c new file mode 100644 index 00000000..99505911 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_moddiv.c @@ -0,0 +1,488 @@ +/* + * Copyright (c) 2018 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* + * In this file, we handle big integers with a custom format, i.e. + * without the usual one-word header. Value is split into 31-bit words, + * each stored in a 32-bit slot (top bit is zero) in little-endian + * order. The length (in words) is provided explicitly. In some cases, + * the value can be negative (using two's complement representation). In + * some cases, the top word is allowed to have a 32th bit. + */ + +/* + * Negate big integer conditionally. The value consists of 'len' words, + * with 31 bits in each word (the top bit of each word should be 0, + * except possibly for the last word). If 'ctl' is 1, the negation is + * computed; otherwise, if 'ctl' is 0, then the value is unchanged. + */ +static void +cond_negate(uint32_t *a, size_t len, uint32_t ctl) +{ + size_t k; + uint32_t cc, xm; + + cc = ctl; + xm = -ctl >> 1; + for (k = 0; k < len; k ++) { + uint32_t aw; + + aw = a[k]; + aw = (aw ^ xm) + cc; + a[k] = aw & 0x7FFFFFFF; + cc = aw >> 31; + } +} + +/* + * Finish modular reduction. Rules on input parameters: + * + * if neg = 1, then -m <= a < 0 + * if neg = 0, then 0 <= a < 2*m + * + * If neg = 0, then the top word of a[] may use 32 bits. + * + * Also, modulus m must be odd. + */ +static void +finish_mod(uint32_t *a, size_t len, const uint32_t *m, uint32_t neg) +{ + size_t k; + uint32_t cc, xm, ym; + + /* + * First pass: compare a (assumed nonnegative) with m. + * Note that if the final word uses the top extra bit, then + * subtracting m must yield a value less than 2^31, since we + * assumed that a < 2*m. + */ + cc = 0; + for (k = 0; k < len; k ++) { + uint32_t aw, mw; + + aw = a[k]; + mw = m[k]; + cc = (aw - mw - cc) >> 31; + } + + /* + * At this point: + * if neg = 1, then we must add m (regardless of cc) + * if neg = 0 and cc = 0, then we must subtract m + * if neg = 0 and cc = 1, then we must do nothing + */ + xm = -neg >> 1; + ym = -(neg | (1 - cc)); + cc = neg; + for (k = 0; k < len; k ++) { + uint32_t aw, mw; + + aw = a[k]; + mw = (m[k] ^ xm) & ym; + aw = aw - mw - cc; + a[k] = aw & 0x7FFFFFFF; + cc = aw >> 31; + } +} + +/* + * Compute: + * a <- (a*pa+b*pb)/(2^31) + * b <- (a*qa+b*qb)/(2^31) + * The division is assumed to be exact (i.e. the low word is dropped). + * If the final a is negative, then it is negated. Similarly for b. + * Returned value is the combination of two bits: + * bit 0: 1 if a had to be negated, 0 otherwise + * bit 1: 1 if b had to be negated, 0 otherwise + * + * Factors pa, pb, qa and qb must be at most 2^31 in absolute value. + * Source integers a and b must be nonnegative; top word is not allowed + * to contain an extra 32th bit. + */ +static uint32_t +co_reduce(uint32_t *a, uint32_t *b, size_t len, + int64_t pa, int64_t pb, int64_t qa, int64_t qb) +{ + size_t k; + int64_t cca, ccb; + uint32_t nega, negb; + + cca = 0; + ccb = 0; + for (k = 0; k < len; k ++) { + uint32_t wa, wb; + uint64_t za, zb; + uint64_t tta, ttb; + + /* + * Since: + * |pa| <= 2^31 + * |pb| <= 2^31 + * 0 <= wa <= 2^31 - 1 + * 0 <= wb <= 2^31 - 1 + * |cca| <= 2^32 - 1 + * Then: + * |za| <= (2^31-1)*(2^32) + (2^32-1) = 2^63 - 1 + * + * Thus, the new value of cca is such that |cca| <= 2^32 - 1. + * The same applies to ccb. + */ + wa = a[k]; + wb = b[k]; + za = wa * (uint64_t)pa + wb * (uint64_t)pb + (uint64_t)cca; + zb = wa * (uint64_t)qa + wb * (uint64_t)qb + (uint64_t)ccb; + if (k > 0) { + a[k - 1] = za & 0x7FFFFFFF; + b[k - 1] = zb & 0x7FFFFFFF; + } + + /* + * For the new values of cca and ccb, we need a signed + * right-shift; since, in C, right-shifting a signed + * negative value is implementation-defined, we use a + * custom portable sign extension expression. + */ +#define M ((uint64_t)1 << 32) + tta = za >> 31; + ttb = zb >> 31; + tta = (tta ^ M) - M; + ttb = (ttb ^ M) - M; + cca = *(int64_t *)&tta; + ccb = *(int64_t *)&ttb; +#undef M + } + a[len - 1] = (uint32_t)cca; + b[len - 1] = (uint32_t)ccb; + + nega = (uint32_t)((uint64_t)cca >> 63); + negb = (uint32_t)((uint64_t)ccb >> 63); + cond_negate(a, len, nega); + cond_negate(b, len, negb); + return nega | (negb << 1); +} + +/* + * Compute: + * a <- (a*pa+b*pb)/(2^31) mod m + * b <- (a*qa+b*qb)/(2^31) mod m + * + * m0i is equal to -1/m[0] mod 2^31. + * + * Factors pa, pb, qa and qb must be at most 2^31 in absolute value. + * Source integers a and b must be nonnegative; top word is not allowed + * to contain an extra 32th bit. + */ +static void +co_reduce_mod(uint32_t *a, uint32_t *b, size_t len, + int64_t pa, int64_t pb, int64_t qa, int64_t qb, + const uint32_t *m, uint32_t m0i) +{ + size_t k; + int64_t cca, ccb; + uint32_t fa, fb; + + cca = 0; + ccb = 0; + fa = ((a[0] * (uint32_t)pa + b[0] * (uint32_t)pb) * m0i) & 0x7FFFFFFF; + fb = ((a[0] * (uint32_t)qa + b[0] * (uint32_t)qb) * m0i) & 0x7FFFFFFF; + for (k = 0; k < len; k ++) { + uint32_t wa, wb; + uint64_t za, zb; + uint64_t tta, ttb; + + /* + * In this loop, carries 'cca' and 'ccb' always fit on + * 33 bits (in absolute value). + */ + wa = a[k]; + wb = b[k]; + za = wa * (uint64_t)pa + wb * (uint64_t)pb + + m[k] * (uint64_t)fa + (uint64_t)cca; + zb = wa * (uint64_t)qa + wb * (uint64_t)qb + + m[k] * (uint64_t)fb + (uint64_t)ccb; + if (k > 0) { + a[k - 1] = (uint32_t)za & 0x7FFFFFFF; + b[k - 1] = (uint32_t)zb & 0x7FFFFFFF; + } + +#define M ((uint64_t)1 << 32) + tta = za >> 31; + ttb = zb >> 31; + tta = (tta ^ M) - M; + ttb = (ttb ^ M) - M; + cca = *(int64_t *)&tta; + ccb = *(int64_t *)&ttb; +#undef M + } + a[len - 1] = (uint32_t)cca; + b[len - 1] = (uint32_t)ccb; + + /* + * At this point: + * -m <= a < 2*m + * -m <= b < 2*m + * (this is a case of Montgomery reduction) + * The top word of 'a' and 'b' may have a 32-th bit set. + * We may have to add or subtract the modulus. + */ + finish_mod(a, len, m, (uint32_t)((uint64_t)cca >> 63)); + finish_mod(b, len, m, (uint32_t)((uint64_t)ccb >> 63)); +} + +/* see inner.h */ +uint32_t +br_i31_moddiv(uint32_t *x, const uint32_t *y, const uint32_t *m, uint32_t m0i, + uint32_t *t) +{ + /* + * Algorithm is an extended binary GCD. We maintain four values + * a, b, u and v, with the following invariants: + * + * a * x = y * u mod m + * b * x = y * v mod m + * + * Starting values are: + * + * a = y + * b = m + * u = x + * v = 0 + * + * The formal definition of the algorithm is a sequence of steps: + * + * - If a is even, then a <- a/2 and u <- u/2 mod m. + * - Otherwise, if b is even, then b <- b/2 and v <- v/2 mod m. + * - Otherwise, if a > b, then a <- (a-b)/2 and u <- (u-v)/2 mod m. + * - Otherwise, b <- (b-a)/2 and v <- (v-u)/2 mod m. + * + * Algorithm stops when a = b. At that point, they both are equal + * to GCD(y,m); the modular division succeeds if that value is 1. + * The result of the modular division is then u (or v: both are + * equal at that point). + * + * Each step makes either a or b shrink by at least one bit; hence, + * if m has bit length k bits, then 2k-2 steps are sufficient. + * + * + * Though complexity is quadratic in the size of m, the bit-by-bit + * processing is not very efficient. We can speed up processing by + * remarking that the decisions are taken based only on observation + * of the top and low bits of a and b. + * + * In the loop below, at each iteration, we use the two top words + * of a and b, and the low words of a and b, to compute reduction + * parameters pa, pb, qa and qb such that the new values for a + * and b are: + * + * a' = (a*pa + b*pb) / (2^31) + * b' = (a*qa + b*qb) / (2^31) + * + * the division being exact. + * + * Since the choices are based on the top words, they may be slightly + * off, requiring an optional correction: if a' < 0, then we replace + * pa with -pa, and pb with -pb. The total length of a and b is + * thus reduced by at least 30 bits at each iteration. + * + * The stopping conditions are still the same, though: when a + * and b become equal, they must be both odd (since m is odd, + * the GCD cannot be even), therefore the next operation is a + * subtraction, and one of the values becomes 0. At that point, + * nothing else happens, i.e. one value is stuck at 0, and the + * other one is the GCD. + */ + size_t len, k; + uint32_t *a, *b, *u, *v; + uint32_t num, r; + + len = (m[0] + 31) >> 5; + a = t; + b = a + len; + u = x + 1; + v = b + len; + memcpy(a, y + 1, len * sizeof *y); + memcpy(b, m + 1, len * sizeof *m); + memset(v, 0, len * sizeof *v); + + /* + * Loop below ensures that a and b are reduced by some bits each, + * for a total of at least 30 bits. + */ + for (num = ((m[0] - (m[0] >> 5)) << 1) + 30; num >= 30; num -= 30) { + size_t j; + uint32_t c0, c1; + uint32_t a0, a1, b0, b1; + uint64_t a_hi, b_hi; + uint32_t a_lo, b_lo; + int64_t pa, pb, qa, qb; + int i; + + /* + * Extract top words of a and b. If j is the highest + * index >= 1 such that a[j] != 0 or b[j] != 0, then we want + * (a[j] << 31) + a[j - 1], and (b[j] << 31) + b[j - 1]. + * If a and b are down to one word each, then we use a[0] + * and b[0]. + */ + c0 = (uint32_t)-1; + c1 = (uint32_t)-1; + a0 = 0; + a1 = 0; + b0 = 0; + b1 = 0; + j = len; + while (j -- > 0) { + uint32_t aw, bw; + + aw = a[j]; + bw = b[j]; + a0 ^= (a0 ^ aw) & c0; + a1 ^= (a1 ^ aw) & c1; + b0 ^= (b0 ^ bw) & c0; + b1 ^= (b1 ^ bw) & c1; + c1 = c0; + c0 &= (((aw | bw) + 0x7FFFFFFF) >> 31) - (uint32_t)1; + } + + /* + * If c1 = 0, then we grabbed two words for a and b. + * If c1 != 0 but c0 = 0, then we grabbed one word. It + * is not possible that c1 != 0 and c0 != 0, because that + * would mean that both integers are zero. + */ + a1 |= a0 & c1; + a0 &= ~c1; + b1 |= b0 & c1; + b0 &= ~c1; + a_hi = ((uint64_t)a0 << 31) + a1; + b_hi = ((uint64_t)b0 << 31) + b1; + a_lo = a[0]; + b_lo = b[0]; + + /* + * Compute reduction factors: + * + * a' = a*pa + b*pb + * b' = a*qa + b*qb + * + * such that a' and b' are both multiple of 2^31, but are + * only marginally larger than a and b. + */ + pa = 1; + pb = 0; + qa = 0; + qb = 1; + for (i = 0; i < 31; i ++) { + /* + * At each iteration: + * + * a <- (a-b)/2 if: a is odd, b is odd, a_hi > b_hi + * b <- (b-a)/2 if: a is odd, b is odd, a_hi <= b_hi + * a <- a/2 if: a is even + * b <- b/2 if: a is odd, b is even + * + * We multiply a_lo and b_lo by 2 at each + * iteration, thus a division by 2 really is a + * non-multiplication by 2. + */ + uint32_t r, oa, ob, cAB, cBA, cA; + uint64_t rz; + + /* + * r = GT(a_hi, b_hi) + * But the GT() function works on uint32_t operands, + * so we inline a 64-bit version here. + */ + rz = b_hi - a_hi; + r = (uint32_t)((rz ^ ((a_hi ^ b_hi) + & (a_hi ^ rz))) >> 63); + + /* + * cAB = 1 if b must be subtracted from a + * cBA = 1 if a must be subtracted from b + * cA = 1 if a is divided by 2, 0 otherwise + * + * Rules: + * + * cAB and cBA cannot be both 1. + * if a is not divided by 2, b is. + */ + oa = (a_lo >> i) & 1; + ob = (b_lo >> i) & 1; + cAB = oa & ob & r; + cBA = oa & ob & NOT(r); + cA = cAB | NOT(oa); + + /* + * Conditional subtractions. + */ + a_lo -= b_lo & -cAB; + a_hi -= b_hi & -(uint64_t)cAB; + pa -= qa & -(int64_t)cAB; + pb -= qb & -(int64_t)cAB; + b_lo -= a_lo & -cBA; + b_hi -= a_hi & -(uint64_t)cBA; + qa -= pa & -(int64_t)cBA; + qb -= pb & -(int64_t)cBA; + + /* + * Shifting. + */ + a_lo += a_lo & (cA - 1); + pa += pa & ((int64_t)cA - 1); + pb += pb & ((int64_t)cA - 1); + a_hi ^= (a_hi ^ (a_hi >> 1)) & -(uint64_t)cA; + b_lo += b_lo & -cA; + qa += qa & -(int64_t)cA; + qb += qb & -(int64_t)cA; + b_hi ^= (b_hi ^ (b_hi >> 1)) & ((uint64_t)cA - 1); + } + + /* + * Replace a and b with new values a' and b'. + */ + r = co_reduce(a, b, len, pa, pb, qa, qb); + pa -= pa * ((r & 1) << 1); + pb -= pb * ((r & 1) << 1); + qa -= qa * (r & 2); + qb -= qb * (r & 2); + co_reduce_mod(u, v, len, pa, pb, qa, qb, m + 1, m0i); + } + + /* + * Now one of the arrays should be 0, and the other contains + * the GCD. If a is 0, then u is 0 as well, and v contains + * the division result. + * Result is correct if and only if GCD is 1. + */ + r = (a[0] | b[0]) ^ 1; + u[0] |= v[0]; + for (k = 1; k < len; k ++) { + r |= a[k] | b[k]; + u[k] |= v[k]; + } + return EQ0(r); +} diff --git a/test/monniaux/BearSSL/src/int/i31_modpow.c b/test/monniaux/BearSSL/src/int/i31_modpow.c new file mode 100644 index 00000000..4ef3f5d5 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_modpow.c @@ -0,0 +1,65 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i31_modpow(uint32_t *x, + const unsigned char *e, size_t elen, + const uint32_t *m, uint32_t m0i, uint32_t *t1, uint32_t *t2) +{ + size_t mlen; + uint32_t k; + + /* + * 'mlen' is the length of m[] expressed in bytes (including + * the "bit length" first field). + */ + mlen = ((m[0] + 63) >> 5) * sizeof m[0]; + + /* + * Throughout the algorithm: + * -- t1[] is in Montgomery representation; it contains x, x^2, + * x^4, x^8... + * -- The result is accumulated, in normal representation, in + * the x[] array. + * -- t2[] is used as destination buffer for each multiplication. + * + * Note that there is no need to call br_i32_from_monty(). + */ + memcpy(t1, x, mlen); + br_i31_to_monty(t1, m); + br_i31_zero(x, m[0]); + x[1] = 1; + for (k = 0; k < ((uint32_t)elen << 3); k ++) { + uint32_t ctl; + + ctl = (e[elen - 1 - (k >> 3)] >> (k & 7)) & 1; + br_i31_montymul(t2, x, t1, m, m0i); + CCOPY(ctl, x, t2, mlen); + br_i31_montymul(t2, t1, t1, m, m0i); + memcpy(t1, t2, mlen); + } +} diff --git a/test/monniaux/BearSSL/src/int/i31_modpow2.c b/test/monniaux/BearSSL/src/int/i31_modpow2.c new file mode 100644 index 00000000..0b8f8cf7 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_modpow2.c @@ -0,0 +1,160 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i31_modpow_opt(uint32_t *x, + const unsigned char *e, size_t elen, + const uint32_t *m, uint32_t m0i, uint32_t *tmp, size_t twlen) +{ + size_t mlen, mwlen; + uint32_t *t1, *t2, *base; + size_t u, v; + uint32_t acc; + int acc_len, win_len; + + /* + * Get modulus size. + */ + mwlen = (m[0] + 63) >> 5; + mlen = mwlen * sizeof m[0]; + mwlen += (mwlen & 1); + t1 = tmp; + t2 = tmp + mwlen; + + /* + * Compute possible window size, with a maximum of 5 bits. + * When the window has size 1 bit, we use a specific code + * that requires only two temporaries. Otherwise, for a + * window of k bits, we need 2^k+1 temporaries. + */ + if (twlen < (mwlen << 1)) { + return 0; + } + for (win_len = 5; win_len > 1; win_len --) { + if ((((uint32_t)1 << win_len) + 1) * mwlen <= twlen) { + break; + } + } + + /* + * Everything is done in Montgomery representation. + */ + br_i31_to_monty(x, m); + + /* + * Compute window contents. If the window has size one bit only, + * then t2 is set to x; otherwise, t2[0] is left untouched, and + * t2[k] is set to x^k (for k >= 1). + */ + if (win_len == 1) { + memcpy(t2, x, mlen); + } else { + memcpy(t2 + mwlen, x, mlen); + base = t2 + mwlen; + for (u = 2; u < ((unsigned)1 << win_len); u ++) { + br_i31_montymul(base + mwlen, base, x, m, m0i); + base += mwlen; + } + } + + /* + * We need to set x to 1, in Montgomery representation. This can + * be done efficiently by setting the high word to 1, then doing + * one word-sized shift. + */ + br_i31_zero(x, m[0]); + x[(m[0] + 31) >> 5] = 1; + br_i31_muladd_small(x, 0, m); + + /* + * We process bits from most to least significant. At each + * loop iteration, we have acc_len bits in acc. + */ + acc = 0; + acc_len = 0; + while (acc_len > 0 || elen > 0) { + int i, k; + uint32_t bits; + + /* + * Get the next bits. + */ + k = win_len; + if (acc_len < win_len) { + if (elen > 0) { + acc = (acc << 8) | *e ++; + elen --; + acc_len += 8; + } else { + k = acc_len; + } + } + bits = (acc >> (acc_len - k)) & (((uint32_t)1 << k) - 1); + acc_len -= k; + + /* + * We could get exactly k bits. Compute k squarings. + */ + for (i = 0; i < k; i ++) { + br_i31_montymul(t1, x, x, m, m0i); + memcpy(x, t1, mlen); + } + + /* + * Window lookup: we want to set t2 to the window + * lookup value, assuming the bits are non-zero. If + * the window length is 1 bit only, then t2 is + * already set; otherwise, we do a constant-time lookup. + */ + if (win_len > 1) { + br_i31_zero(t2, m[0]); + base = t2 + mwlen; + for (u = 1; u < ((uint32_t)1 << k); u ++) { + uint32_t mask; + + mask = -EQ(u, bits); + for (v = 1; v < mwlen; v ++) { + t2[v] |= mask & base[v]; + } + base += mwlen; + } + } + + /* + * Multiply with the looked-up value. We keep the + * product only if the exponent bits are not all-zero. + */ + br_i31_montymul(t1, x, t2, m, m0i); + CCOPY(NEQ(bits, 0), x, t1, mlen); + } + + /* + * Convert back from Montgomery representation, and exit. + */ + br_i31_from_monty(x, m, m0i); + return 1; +} diff --git a/test/monniaux/BearSSL/src/int/i31_montmul.c b/test/monniaux/BearSSL/src/int/i31_montmul.c new file mode 100644 index 00000000..758f8f42 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_montmul.c @@ -0,0 +1,127 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i31_montymul(uint32_t *d, const uint32_t *x, const uint32_t *y, + const uint32_t *m, uint32_t m0i) +{ + /* + * Each outer loop iteration computes: + * d <- (d + xu*y + f*m) / 2^31 + * We have xu <= 2^31-1 and f <= 2^31-1. + * Thus, if d <= 2*m-1 on input, then: + * 2*m-1 + 2*(2^31-1)*m <= (2^32)*m-1 + * and the new d value is less than 2*m. + * + * We represent d over 31-bit words, with an extra word 'dh' + * which can thus be only 0 or 1. + */ + size_t len, len4, u, v; + uint32_t dh; + + len = (m[0] + 31) >> 5; + len4 = len & ~(size_t)3; + br_i31_zero(d, m[0]); + dh = 0; + for (u = 0; u < len; u ++) { + /* + * The carry for each operation fits on 32 bits: + * d[v+1] <= 2^31-1 + * xu*y[v+1] <= (2^31-1)*(2^31-1) + * f*m[v+1] <= (2^31-1)*(2^31-1) + * r <= 2^32-1 + * (2^31-1) + 2*(2^31-1)*(2^31-1) + (2^32-1) = 2^63 - 2^31 + * After division by 2^31, the new r is then at most 2^32-1 + * + * Using a 32-bit carry has performance benefits on 32-bit + * systems; however, on 64-bit architectures, we prefer to + * keep the carry (r) in a 64-bit register, thus avoiding some + * "clear high bits" operations. + */ + uint32_t f, xu; +#if BR_64 + uint64_t r; +#else + uint32_t r; +#endif + + xu = x[u + 1]; + f = MUL31_lo((d[1] + MUL31_lo(x[u + 1], y[1])), m0i); + + r = 0; + for (v = 0; v < len4; v += 4) { + uint64_t z; + + z = (uint64_t)d[v + 1] + MUL31(xu, y[v + 1]) + + MUL31(f, m[v + 1]) + r; + r = z >> 31; + d[v + 0] = (uint32_t)z & 0x7FFFFFFF; + z = (uint64_t)d[v + 2] + MUL31(xu, y[v + 2]) + + MUL31(f, m[v + 2]) + r; + r = z >> 31; + d[v + 1] = (uint32_t)z & 0x7FFFFFFF; + z = (uint64_t)d[v + 3] + MUL31(xu, y[v + 3]) + + MUL31(f, m[v + 3]) + r; + r = z >> 31; + d[v + 2] = (uint32_t)z & 0x7FFFFFFF; + z = (uint64_t)d[v + 4] + MUL31(xu, y[v + 4]) + + MUL31(f, m[v + 4]) + r; + r = z >> 31; + d[v + 3] = (uint32_t)z & 0x7FFFFFFF; + } + for (; v < len; v ++) { + uint64_t z; + + z = (uint64_t)d[v + 1] + MUL31(xu, y[v + 1]) + + MUL31(f, m[v + 1]) + r; + r = z >> 31; + d[v] = (uint32_t)z & 0x7FFFFFFF; + } + + /* + * Since the new dh can only be 0 or 1, the addition of + * the old dh with the carry MUST fit on 32 bits, and + * thus can be done into dh itself. + */ + dh += r; + d[len] = dh & 0x7FFFFFFF; + dh >>= 31; + } + + /* + * We must write back the bit length because it was overwritten in + * the loop (not overwriting it would require a test in the loop, + * which would yield bigger and slower code). + */ + d[0] = m[0]; + + /* + * d[] may still be greater than m[] at that point; notably, the + * 'dh' word may be non-zero. + */ + br_i31_sub(d, m, NEQ(dh, 0) | NOT(br_i31_sub(d, m, 0))); +} diff --git a/test/monniaux/BearSSL/src/int/i31_mulacc.c b/test/monniaux/BearSSL/src/int/i31_mulacc.c new file mode 100644 index 00000000..7410e546 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_mulacc.c @@ -0,0 +1,74 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i31_mulacc(uint32_t *d, const uint32_t *a, const uint32_t *b) +{ + size_t alen, blen, u; + uint32_t dl, dh; + + alen = (a[0] + 31) >> 5; + blen = (b[0] + 31) >> 5; + + /* + * We want to add the two bit lengths, but these are encoded, + * which requires some extra care. + */ + dl = (a[0] & 31) + (b[0] & 31); + dh = (a[0] >> 5) + (b[0] >> 5); + d[0] = (dh << 5) + dl + (~(uint32_t)(dl - 31) >> 31); + + for (u = 0; u < blen; u ++) { + uint32_t f; + size_t v; + + /* + * Carry always fits on 31 bits; we want to keep it in a + * 32-bit register on 32-bit architectures (on a 64-bit + * architecture, cast down from 64 to 32 bits means + * clearing the high bits, which is not free; on a 32-bit + * architecture, the same operation really means ignoring + * the top register, which has negative or zero cost). + */ +#if BR_64 + uint64_t cc; +#else + uint32_t cc; +#endif + + f = b[1 + u]; + cc = 0; + for (v = 0; v < alen; v ++) { + uint64_t z; + + z = (uint64_t)d[1 + u + v] + MUL31(f, a[1 + v]) + cc; + cc = z >> 31; + d[1 + u + v] = (uint32_t)z & 0x7FFFFFFF; + } + d[1 + u + alen] = (uint32_t)cc; + } +} diff --git a/test/monniaux/BearSSL/src/int/i31_muladd.c b/test/monniaux/BearSSL/src/int/i31_muladd.c new file mode 100644 index 00000000..eecd9e2c --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_muladd.c @@ -0,0 +1,157 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i31_muladd_small(uint32_t *x, uint32_t z, const uint32_t *m) +{ + uint32_t m_bitlen; + unsigned mblr; + size_t u, mlen; + uint32_t a0, a1, b0, hi, g, q, tb; + uint32_t under, over; + uint32_t cc; + + /* + * We can test on the modulus bit length since we accept to + * leak that length. + */ + m_bitlen = m[0]; + if (m_bitlen == 0) { + return; + } + if (m_bitlen <= 31) { + uint32_t lo; + + hi = x[1] >> 1; + lo = (x[1] << 31) | z; + x[1] = br_rem(hi, lo, m[1]); + return; + } + mlen = (m_bitlen + 31) >> 5; + mblr = (unsigned)m_bitlen & 31; + + /* + * Principle: we estimate the quotient (x*2^31+z)/m by + * doing a 64/32 division with the high words. + * + * Let: + * w = 2^31 + * a = (w*a0 + a1) * w^N + a2 + * b = b0 * w^N + b2 + * such that: + * 0 <= a0 < w + * 0 <= a1 < w + * 0 <= a2 < w^N + * w/2 <= b0 < w + * 0 <= b2 < w^N + * a < w*b + * I.e. the two top words of a are a0:a1, the top word of b is + * b0, we ensured that b0 is "full" (high bit set), and a is + * such that the quotient q = a/b fits on one word (0 <= q < w). + * + * If a = b*q + r (with 0 <= r < q), we can estimate q by + * doing an Euclidean division on the top words: + * a0*w+a1 = b0*u + v (with 0 <= v < b0) + * Then the following holds: + * 0 <= u <= w + * u-2 <= q <= u + */ + hi = x[mlen]; + if (mblr == 0) { + a0 = x[mlen]; + memmove(x + 2, x + 1, (mlen - 1) * sizeof *x); + x[1] = z; + a1 = x[mlen]; + b0 = m[mlen]; + } else { + a0 = ((x[mlen] << (31 - mblr)) | (x[mlen - 1] >> mblr)) + & 0x7FFFFFFF; + memmove(x + 2, x + 1, (mlen - 1) * sizeof *x); + x[1] = z; + a1 = ((x[mlen] << (31 - mblr)) | (x[mlen - 1] >> mblr)) + & 0x7FFFFFFF; + b0 = ((m[mlen] << (31 - mblr)) | (m[mlen - 1] >> mblr)) + & 0x7FFFFFFF; + } + + /* + * We estimate a divisor q. If the quotient returned by br_div() + * is g: + * -- If a0 == b0 then g == 0; we want q = 0x7FFFFFFF. + * -- Otherwise: + * -- if g == 0 then we set q = 0; + * -- otherwise, we set q = g - 1. + * The properties described above then ensure that the true + * quotient is q-1, q or q+1. + * + * Take care that a0, a1 and b0 are 31-bit words, not 32-bit. We + * must adjust the parameters to br_div() accordingly. + */ + g = br_div(a0 >> 1, a1 | (a0 << 31), b0); + q = MUX(EQ(a0, b0), 0x7FFFFFFF, MUX(EQ(g, 0), 0, g - 1)); + + /* + * We subtract q*m from x (with the extra high word of value 'hi'). + * Since q may be off by 1 (in either direction), we may have to + * add or subtract m afterwards. + * + * The 'tb' flag will be true (1) at the end of the loop if the + * result is greater than or equal to the modulus (not counting + * 'hi' or the carry). + */ + cc = 0; + tb = 1; + for (u = 1; u <= mlen; u ++) { + uint32_t mw, zw, xw, nxw; + uint64_t zl; + + mw = m[u]; + zl = MUL31(mw, q) + cc; + cc = (uint32_t)(zl >> 31); + zw = (uint32_t)zl & (uint32_t)0x7FFFFFFF; + xw = x[u]; + nxw = xw - zw; + cc += nxw >> 31; + nxw &= 0x7FFFFFFF; + x[u] = nxw; + tb = MUX(EQ(nxw, mw), tb, GT(nxw, mw)); + } + + /* + * If we underestimated q, then either cc < hi (one extra bit + * beyond the top array word), or cc == hi and tb is true (no + * extra bit, but the result is not lower than the modulus). In + * these cases we must subtract m once. + * + * Otherwise, we may have overestimated, which will show as + * cc > hi (thus a negative result). Correction is adding m once. + */ + over = GT(cc, hi); + under = ~over & (tb | LT(cc, hi)); + br_i31_add(x, m, over); + br_i31_sub(x, m, under); +} diff --git a/test/monniaux/BearSSL/src/int/i31_ninv31.c b/test/monniaux/BearSSL/src/int/i31_ninv31.c new file mode 100644 index 00000000..dd83c96a --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_ninv31.c @@ -0,0 +1,39 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i31_ninv31(uint32_t x) +{ + uint32_t y; + + y = 2 - x; + y *= 2 - y * x; + y *= 2 - y * x; + y *= 2 - y * x; + y *= 2 - y * x; + return MUX(x & 1, -y, 0) & 0x7FFFFFFF; +} diff --git a/test/monniaux/BearSSL/src/int/i31_reduce.c b/test/monniaux/BearSSL/src/int/i31_reduce.c new file mode 100644 index 00000000..5c9523ed --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_reduce.c @@ -0,0 +1,66 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i31_reduce(uint32_t *x, const uint32_t *a, const uint32_t *m) +{ + uint32_t m_bitlen, a_bitlen; + size_t mlen, alen, u; + + m_bitlen = m[0]; + mlen = (m_bitlen + 31) >> 5; + + x[0] = m_bitlen; + if (m_bitlen == 0) { + return; + } + + /* + * If the source is shorter, then simply copy all words from a[] + * and zero out the upper words. + */ + a_bitlen = a[0]; + alen = (a_bitlen + 31) >> 5; + if (a_bitlen < m_bitlen) { + memcpy(x + 1, a + 1, alen * sizeof *a); + for (u = alen; u < mlen; u ++) { + x[u + 1] = 0; + } + return; + } + + /* + * The source length is at least equal to that of the modulus. + * We must thus copy N-1 words, and input the remaining words + * one by one. + */ + memcpy(x + 1, a + 2 + (alen - mlen), (mlen - 1) * sizeof *a); + x[mlen] = 0; + for (u = 1 + alen - mlen; u > 0; u --) { + br_i31_muladd_small(x, a[u], m); + } +} diff --git a/test/monniaux/BearSSL/src/int/i31_rshift.c b/test/monniaux/BearSSL/src/int/i31_rshift.c new file mode 100644 index 00000000..db6ba0b8 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_rshift.c @@ -0,0 +1,47 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i31_rshift(uint32_t *x, int count) +{ + size_t u, len; + uint32_t r; + + len = (x[0] + 31) >> 5; + if (len == 0) { + return; + } + r = x[1] >> count; + for (u = 2; u <= len; u ++) { + uint32_t w; + + w = x[u]; + x[u - 1] = ((w << (31 - count)) | r) & 0x7FFFFFFF; + r = w >> count; + } + x[len] = r; +} diff --git a/test/monniaux/BearSSL/src/int/i31_sub.c b/test/monniaux/BearSSL/src/int/i31_sub.c new file mode 100644 index 00000000..39108951 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_sub.c @@ -0,0 +1,46 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i31_sub(uint32_t *a, const uint32_t *b, uint32_t ctl) +{ + uint32_t cc; + size_t u, m; + + cc = 0; + m = (a[0] + 63) >> 5; + for (u = 1; u < m; u ++) { + uint32_t aw, bw, naw; + + aw = a[u]; + bw = b[u]; + naw = aw - bw - cc; + cc = naw >> 31; + a[u] = MUX(ctl, naw & 0x7FFFFFFF, aw); + } + return cc; +} diff --git a/test/monniaux/BearSSL/src/int/i31_tmont.c b/test/monniaux/BearSSL/src/int/i31_tmont.c new file mode 100644 index 00000000..4798ff65 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i31_tmont.c @@ -0,0 +1,36 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i31_to_monty(uint32_t *x, const uint32_t *m) +{ + uint32_t k; + + for (k = (m[0] + 31) >> 5; k > 0; k --) { + br_i31_muladd_small(x, 0, m); + } +} diff --git a/test/monniaux/BearSSL/src/int/i32_add.c b/test/monniaux/BearSSL/src/int/i32_add.c new file mode 100644 index 00000000..620baffd --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_add.c @@ -0,0 +1,51 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i32_add(uint32_t *a, const uint32_t *b, uint32_t ctl) +{ + uint32_t cc; + size_t u, m; + + cc = 0; + m = (a[0] + 63) >> 5; + for (u = 1; u < m; u ++) { + uint32_t aw, bw, naw; + + aw = a[u]; + bw = b[u]; + naw = aw + bw + cc; + + /* + * Carry is 1 if naw < aw. Carry is also 1 if naw == aw + * AND the carry was already 1. + */ + cc = (cc & EQ(naw, aw)) | LT(naw, aw); + a[u] = MUX(ctl, naw, aw); + } + return cc; +} diff --git a/test/monniaux/BearSSL/src/int/i32_bitlen.c b/test/monniaux/BearSSL/src/int/i32_bitlen.c new file mode 100644 index 00000000..40ce9fa0 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_bitlen.c @@ -0,0 +1,44 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i32_bit_length(uint32_t *x, size_t xlen) +{ + uint32_t tw, twk; + + tw = 0; + twk = 0; + while (xlen -- > 0) { + uint32_t w, c; + + c = EQ(tw, 0); + w = x[xlen]; + tw = MUX(c, w, tw); + twk = MUX(c, (uint32_t)xlen, twk); + } + return (twk << 5) + BIT_LENGTH(tw); +} diff --git a/test/monniaux/BearSSL/src/int/i32_decmod.c b/test/monniaux/BearSSL/src/int/i32_decmod.c new file mode 100644 index 00000000..a859af12 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_decmod.c @@ -0,0 +1,77 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i32_decode_mod(uint32_t *x, const void *src, size_t len, const uint32_t *m) +{ + const unsigned char *buf; + uint32_t r; + size_t u, v, mlen; + + buf = src; + + /* + * First pass: determine whether the value fits. The 'r' value + * will contain the comparison result, as 0x00000000 (value is + * equal to the modulus), 0x00000001 (value is greater than the + * modulus), or 0xFFFFFFFF (value is lower than the modulus). + */ + mlen = (m[0] + 7) >> 3; + r = 0; + for (u = (mlen > len) ? mlen : len; u > 0; u --) { + uint32_t mb, xb; + + v = u - 1; + if (v >= mlen) { + mb = 0; + } else { + mb = (m[1 + (v >> 2)] >> ((v & 3) << 3)) & 0xFF; + } + if (v >= len) { + xb = 0; + } else { + xb = buf[len - u]; + } + r = MUX(EQ(r, 0), (uint32_t)CMP(xb, mb), r); + } + + /* + * Only r == 0xFFFFFFFF is acceptable. We want to set r to 0xFF if + * the value fits, 0x00 otherwise. + */ + r >>= 24; + br_i32_zero(x, m[0]); + u = (mlen > len) ? len : mlen; + while (u > 0) { + uint32_t xb; + + xb = buf[len - u] & r; + u --; + x[1 + (u >> 2)] |= xb << ((u & 3) << 3); + } + return r >> 7; +} diff --git a/test/monniaux/BearSSL/src/int/i32_decode.c b/test/monniaux/BearSSL/src/int/i32_decode.c new file mode 100644 index 00000000..f2890384 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_decode.c @@ -0,0 +1,63 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i32_decode(uint32_t *x, const void *src, size_t len) +{ + const unsigned char *buf; + size_t u, v; + + buf = src; + u = len; + v = 1; + for (;;) { + if (u < 4) { + uint32_t w; + + if (u < 2) { + if (u == 0) { + break; + } else { + w = buf[0]; + } + } else { + if (u == 2) { + w = br_dec16be(buf); + } else { + w = ((uint32_t)buf[0] << 16) + | br_dec16be(buf + 1); + } + } + x[v ++] = w; + break; + } else { + u -= 4; + x[v ++] = br_dec32be(buf + u); + } + } + x[0] = br_i32_bit_length(x + 1, v - 1); +} diff --git a/test/monniaux/BearSSL/src/int/i32_decred.c b/test/monniaux/BearSSL/src/int/i32_decred.c new file mode 100644 index 00000000..dc476db0 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_decred.c @@ -0,0 +1,107 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i32_decode_reduce(uint32_t *x, + const void *src, size_t len, const uint32_t *m) +{ + uint32_t m_bitlen; + size_t mblen, k, q; + const unsigned char *buf; + + m_bitlen = m[0]; + + /* + * Special case for an invalid modulus. + */ + if (m_bitlen == 0) { + x[0] = 0; + return; + } + + /* + * Clear the destination. + */ + br_i32_zero(x, m_bitlen); + + /* + * First decode directly as many bytes as possible without + * reduction, taking care to leave a number of bytes which + * is a multiple of 4. + */ + mblen = (m_bitlen + 7) >> 3; + k = mblen - 1; + + /* + * Up to k bytes can be safely decoded. + */ + if (k >= len) { + br_i32_decode(x, src, len); + x[0] = m_bitlen; + return; + } + + /* + * We want to first inject some bytes with direct decoding, + * then extra bytes by whole 32-bit words. First compute + * the size that should be injected that way. + */ + buf = src; + q = (len - k + 3) & ~(size_t)3; + + /* + * It may happen that this is more than what we already have + * (by at most 3 bytes). Such a case may happen only with + * a very short modulus. In that case, we must process the first + * bytes "manually". + */ + if (q > len) { + int i; + uint32_t w; + + w = 0; + for (i = 0; i < 4; i ++) { + w <<= 8; + if (q <= len) { + w |= buf[len - q]; + } + q --; + } + br_i32_muladd_small(x, w, m); + } else { + br_i32_decode(x, buf, len - q); + x[0] = m_bitlen; + } + + /* + * At that point, we have exactly q bytes to inject, and q is + * a multiple of 4. + */ + for (k = len - q; k < len; k += 4) { + br_i32_muladd_small(x, br_dec32be(buf + k), m); + } +} diff --git a/test/monniaux/BearSSL/src/int/i32_div32.c b/test/monniaux/BearSSL/src/int/i32_div32.c new file mode 100644 index 00000000..d8b8023d --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_div32.c @@ -0,0 +1,56 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_divrem(uint32_t hi, uint32_t lo, uint32_t d, uint32_t *r) +{ + /* TODO: optimize this */ + uint32_t q; + uint32_t ch, cf; + int k; + + q = 0; + ch = EQ(hi, d); + hi = MUX(ch, 0, hi); + for (k = 31; k > 0; k --) { + int j; + uint32_t w, ctl, hi2, lo2; + + j = 32 - k; + w = (hi << j) | (lo >> k); + ctl = GE(w, d) | (hi >> k); + hi2 = (w - d) >> j; + lo2 = lo - (d << k); + hi = MUX(ctl, hi2, hi); + lo = MUX(ctl, lo2, lo); + q |= ctl << k; + } + cf = GE(lo, d) | hi; + q |= cf; + *r = MUX(cf, lo - d, lo); + return q; +} diff --git a/test/monniaux/BearSSL/src/int/i32_encode.c b/test/monniaux/BearSSL/src/int/i32_encode.c new file mode 100644 index 00000000..303652f9 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_encode.c @@ -0,0 +1,72 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i32_encode(void *dst, size_t len, const uint32_t *x) +{ + unsigned char *buf; + size_t k; + + buf = dst; + + /* + * Compute the announced size of x in bytes; extra bytes are + * filled with zeros. + */ + k = (x[0] + 7) >> 3; + while (len > k) { + *buf ++ = 0; + len --; + } + + /* + * Now we use k as index within x[]. That index starts at 1; + * we initialize it to the topmost complete word, and process + * any remaining incomplete word. + */ + k = (len + 3) >> 2; + switch (len & 3) { + case 3: + *buf ++ = x[k] >> 16; + /* fall through */ + case 2: + *buf ++ = x[k] >> 8; + /* fall through */ + case 1: + *buf ++ = x[k]; + k --; + } + + /* + * Encode all complete words. + */ + while (k > 0) { + br_enc32be(buf, x[k]); + k --; + buf += 4; + } +} diff --git a/test/monniaux/BearSSL/src/int/i32_fmont.c b/test/monniaux/BearSSL/src/int/i32_fmont.c new file mode 100644 index 00000000..dc1c9344 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_fmont.c @@ -0,0 +1,60 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i32_from_monty(uint32_t *x, const uint32_t *m, uint32_t m0i) +{ + size_t len, u, v; + + len = (m[0] + 31) >> 5; + for (u = 0; u < len; u ++) { + uint32_t f; + uint64_t cc; + + f = x[1] * m0i; + cc = 0; + for (v = 0; v < len; v ++) { + uint64_t z; + + z = (uint64_t)x[v + 1] + MUL(f, m[v + 1]) + cc; + cc = z >> 32; + if (v != 0) { + x[v] = (uint32_t)z; + } + } + x[len] = (uint32_t)cc; + } + + /* + * We may have to do an extra subtraction, but only if the + * value in x[] is indeed greater than or equal to that of m[], + * which is why we must do two calls (first call computes the + * carry, second call performs the subtraction only if the carry + * is 0). + */ + br_i32_sub(x, m, NOT(br_i32_sub(x, m, 0))); +} diff --git a/test/monniaux/BearSSL/src/int/i32_iszero.c b/test/monniaux/BearSSL/src/int/i32_iszero.c new file mode 100644 index 00000000..659df7f2 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_iszero.c @@ -0,0 +1,39 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i32_iszero(const uint32_t *x) +{ + uint32_t z; + size_t u; + + z = 0; + for (u = (x[0] + 31) >> 5; u > 0; u --) { + z |= x[u]; + } + return ~(z | -z) >> 31; +} diff --git a/test/monniaux/BearSSL/src/int/i32_modpow.c b/test/monniaux/BearSSL/src/int/i32_modpow.c new file mode 100644 index 00000000..034aba06 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_modpow.c @@ -0,0 +1,65 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i32_modpow(uint32_t *x, + const unsigned char *e, size_t elen, + const uint32_t *m, uint32_t m0i, uint32_t *t1, uint32_t *t2) +{ + size_t mlen; + uint32_t k; + + /* + * 'mlen' is the length of m[] expressed in bytes (including + * the "bit length" first field). + */ + mlen = ((m[0] + 63) >> 5) * sizeof m[0]; + + /* + * Throughout the algorithm: + * -- t1[] is in Montgomery representation; it contains x, x^2, + * x^4, x^8... + * -- The result is accumulated, in normal representation, in + * the x[] array. + * -- t2[] is used as destination buffer for each multiplication. + * + * Note that there is no need to call br_i32_from_monty(). + */ + memcpy(t1, x, mlen); + br_i32_to_monty(t1, m); + br_i32_zero(x, m[0]); + x[1] = 1; + for (k = 0; k < ((uint32_t)elen << 3); k ++) { + uint32_t ctl; + + ctl = (e[elen - 1 - (k >> 3)] >> (k & 7)) & 1; + br_i32_montymul(t2, x, t1, m, m0i); + CCOPY(ctl, x, t2, mlen); + br_i32_montymul(t2, t1, t1, m, m0i); + memcpy(t1, t2, mlen); + } +} diff --git a/test/monniaux/BearSSL/src/int/i32_montmul.c b/test/monniaux/BearSSL/src/int/i32_montmul.c new file mode 100644 index 00000000..7edb376c --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_montmul.c @@ -0,0 +1,69 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i32_montymul(uint32_t *d, const uint32_t *x, const uint32_t *y, + const uint32_t *m, uint32_t m0i) +{ + size_t len, u, v; + uint64_t dh; + + len = (m[0] + 31) >> 5; + br_i32_zero(d, m[0]); + dh = 0; + for (u = 0; u < len; u ++) { + uint32_t f, xu; + uint64_t r1, r2, zh; + + xu = x[u + 1]; + f = (d[1] + x[u + 1] * y[1]) * m0i; + r1 = 0; + r2 = 0; + for (v = 0; v < len; v ++) { + uint64_t z; + uint32_t t; + + z = (uint64_t)d[v + 1] + MUL(xu, y[v + 1]) + r1; + r1 = z >> 32; + t = (uint32_t)z; + z = (uint64_t)t + MUL(f, m[v + 1]) + r2; + r2 = z >> 32; + if (v != 0) { + d[v] = (uint32_t)z; + } + } + zh = dh + r1 + r2; + d[len] = (uint32_t)zh; + dh = zh >> 32; + } + + /* + * d[] may still be greater than m[] at that point; notably, the + * 'dh' word may be non-zero. + */ + br_i32_sub(d, m, NEQ(dh, 0) | NOT(br_i32_sub(d, m, 0))); +} diff --git a/test/monniaux/BearSSL/src/int/i32_mulacc.c b/test/monniaux/BearSSL/src/int/i32_mulacc.c new file mode 100644 index 00000000..55da3858 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_mulacc.c @@ -0,0 +1,56 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i32_mulacc(uint32_t *d, const uint32_t *a, const uint32_t *b) +{ + size_t alen, blen, u; + + alen = (a[0] + 31) >> 5; + blen = (b[0] + 31) >> 5; + d[0] = a[0] + b[0]; + for (u = 0; u < blen; u ++) { + uint32_t f; + size_t v; +#if BR_64 + uint64_t cc; +#else + uint32_t cc; +#endif + + f = b[1 + u]; + cc = 0; + for (v = 0; v < alen; v ++) { + uint64_t z; + + z = (uint64_t)d[1 + u + v] + MUL(f, a[1 + v]) + cc; + cc = z >> 32; + d[1 + u + v] = (uint32_t)z; + } + d[1 + u + alen] = (uint32_t)cc; + } +} diff --git a/test/monniaux/BearSSL/src/int/i32_muladd.c b/test/monniaux/BearSSL/src/int/i32_muladd.c new file mode 100644 index 00000000..dd526ad5 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_muladd.c @@ -0,0 +1,138 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i32_muladd_small(uint32_t *x, uint32_t z, const uint32_t *m) +{ + uint32_t m_bitlen; + size_t u, mlen; + uint32_t a0, a1, b0, hi, g, q, tb; + uint32_t chf, clow, under, over; + uint64_t cc; + + /* + * We can test on the modulus bit length since we accept to + * leak that length. + */ + m_bitlen = m[0]; + if (m_bitlen == 0) { + return; + } + if (m_bitlen <= 32) { + x[1] = br_rem(x[1], z, m[1]); + return; + } + mlen = (m_bitlen + 31) >> 5; + + /* + * Principle: we estimate the quotient (x*2^32+z)/m by + * doing a 64/32 division with the high words. + * + * Let: + * w = 2^32 + * a = (w*a0 + a1) * w^N + a2 + * b = b0 * w^N + b2 + * such that: + * 0 <= a0 < w + * 0 <= a1 < w + * 0 <= a2 < w^N + * w/2 <= b0 < w + * 0 <= b2 < w^N + * a < w*b + * I.e. the two top words of a are a0:a1, the top word of b is + * b0, we ensured that b0 is "full" (high bit set), and a is + * such that the quotient q = a/b fits on one word (0 <= q < w). + * + * If a = b*q + r (with 0 <= r < q), we can estimate q by + * doing an Euclidean division on the top words: + * a0*w+a1 = b0*u + v (with 0 <= v < w) + * Then the following holds: + * 0 <= u <= w + * u-2 <= q <= u + */ + a0 = br_i32_word(x, m_bitlen - 32); + hi = x[mlen]; + memmove(x + 2, x + 1, (mlen - 1) * sizeof *x); + x[1] = z; + a1 = br_i32_word(x, m_bitlen - 32); + b0 = br_i32_word(m, m_bitlen - 32); + + /* + * We estimate a divisor q. If the quotient returned by br_div() + * is g: + * -- If a0 == b0 then g == 0; we want q = 0xFFFFFFFF. + * -- Otherwise: + * -- if g == 0 then we set q = 0; + * -- otherwise, we set q = g - 1. + * The properties described above then ensure that the true + * quotient is q-1, q or q+1. + */ + g = br_div(a0, a1, b0); + q = MUX(EQ(a0, b0), 0xFFFFFFFF, MUX(EQ(g, 0), 0, g - 1)); + + /* + * We subtract q*m from x (with the extra high word of value 'hi'). + * Since q may be off by 1 (in either direction), we may have to + * add or subtract m afterwards. + * + * The 'tb' flag will be true (1) at the end of the loop if the + * result is greater than or equal to the modulus (not counting + * 'hi' or the carry). + */ + cc = 0; + tb = 1; + for (u = 1; u <= mlen; u ++) { + uint32_t mw, zw, xw, nxw; + uint64_t zl; + + mw = m[u]; + zl = MUL(mw, q) + cc; + cc = (uint32_t)(zl >> 32); + zw = (uint32_t)zl; + xw = x[u]; + nxw = xw - zw; + cc += (uint64_t)GT(nxw, xw); + x[u] = nxw; + tb = MUX(EQ(nxw, mw), tb, GT(nxw, mw)); + } + + /* + * If we underestimated q, then either cc < hi (one extra bit + * beyond the top array word), or cc == hi and tb is true (no + * extra bit, but the result is not lower than the modulus). In + * these cases we must subtract m once. + * + * Otherwise, we may have overestimated, which will show as + * cc > hi (thus a negative result). Correction is adding m once. + */ + chf = (uint32_t)(cc >> 32); + clow = (uint32_t)cc; + over = chf | GT(clow, hi); + under = ~over & (tb | (~chf & LT(clow, hi))); + br_i32_add(x, m, over); + br_i32_sub(x, m, under); +} diff --git a/test/monniaux/BearSSL/src/int/i32_ninv32.c b/test/monniaux/BearSSL/src/int/i32_ninv32.c new file mode 100644 index 00000000..65644341 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_ninv32.c @@ -0,0 +1,39 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i32_ninv32(uint32_t x) +{ + uint32_t y; + + y = 2 - x; + y *= 2 - y * x; + y *= 2 - y * x; + y *= 2 - y * x; + y *= 2 - y * x; + return MUX(x & 1, -y, 0); +} diff --git a/test/monniaux/BearSSL/src/int/i32_reduce.c b/test/monniaux/BearSSL/src/int/i32_reduce.c new file mode 100644 index 00000000..90fff092 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_reduce.c @@ -0,0 +1,66 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i32_reduce(uint32_t *x, const uint32_t *a, const uint32_t *m) +{ + uint32_t m_bitlen, a_bitlen; + size_t mlen, alen, u; + + m_bitlen = m[0]; + mlen = (m_bitlen + 31) >> 5; + + x[0] = m_bitlen; + if (m_bitlen == 0) { + return; + } + + /* + * If the source is shorter, then simply copy all words from a[] + * and zero out the upper words. + */ + a_bitlen = a[0]; + alen = (a_bitlen + 31) >> 5; + if (a_bitlen < m_bitlen) { + memcpy(x + 1, a + 1, alen * sizeof *a); + for (u = alen; u < mlen; u ++) { + x[u + 1] = 0; + } + return; + } + + /* + * The source length is at least equal to that of the modulus. + * We must thus copy N-1 words, and input the remaining words + * one by one. + */ + memcpy(x + 1, a + 2 + (alen - mlen), (mlen - 1) * sizeof *a); + x[mlen] = 0; + for (u = 1 + alen - mlen; u > 0; u --) { + br_i32_muladd_small(x, a[u], m); + } +} diff --git a/test/monniaux/BearSSL/src/int/i32_sub.c b/test/monniaux/BearSSL/src/int/i32_sub.c new file mode 100644 index 00000000..9c500238 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_sub.c @@ -0,0 +1,51 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +uint32_t +br_i32_sub(uint32_t *a, const uint32_t *b, uint32_t ctl) +{ + uint32_t cc; + size_t u, m; + + cc = 0; + m = (a[0] + 63) >> 5; + for (u = 1; u < m; u ++) { + uint32_t aw, bw, naw; + + aw = a[u]; + bw = b[u]; + naw = aw - bw - cc; + + /* + * Carry is 1 if naw > aw. Carry is 1 also if naw == aw + * AND the carry was already 1. + */ + cc = (cc & EQ(naw, aw)) | GT(naw, aw); + a[u] = MUX(ctl, naw, aw); + } + return cc; +} diff --git a/test/monniaux/BearSSL/src/int/i32_tmont.c b/test/monniaux/BearSSL/src/int/i32_tmont.c new file mode 100644 index 00000000..058cd886 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i32_tmont.c @@ -0,0 +1,36 @@ +/* + * Copyright (c) 2016 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +/* see inner.h */ +void +br_i32_to_monty(uint32_t *x, const uint32_t *m) +{ + uint32_t k; + + for (k = (m[0] + 31) >> 5; k > 0; k --) { + br_i32_muladd_small(x, 0, m); + } +} diff --git a/test/monniaux/BearSSL/src/int/i62_modpow2.c b/test/monniaux/BearSSL/src/int/i62_modpow2.c new file mode 100644 index 00000000..2db537f0 --- /dev/null +++ b/test/monniaux/BearSSL/src/int/i62_modpow2.c @@ -0,0 +1,493 @@ +/* + * Copyright (c) 2017 Thomas Pornin <pornin@bolet.org> + * + * Permission is hereby granted, free of charge, to any person obtaining + * a copy of this software and associated documentation files (the + * "Software"), to deal in the Software without restriction, including + * without limitation the rights to use, copy, modify, merge, publish, + * distribute, sublicense, and/or sell copies of the Software, and to + * permit persons to whom the Software is furnished to do so, subject to + * the following conditions: + * + * The above copyright notice and this permission notice shall be + * included in all copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, + * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF + * MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND + * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS + * BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN + * ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN + * CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ + +#include "inner.h" + +#if BR_INT128 || BR_UMUL128 + +#if BR_INT128 + +/* + * Compute x*y+v1+v2. Operands are 64-bit, and result is 128-bit, with + * high word in "hi" and low word in "lo". + */ +#define FMA1(hi, lo, x, y, v1, v2) do { \ + unsigned __int128 fmaz; \ + fmaz = (unsigned __int128)(x) * (unsigned __int128)(y) \ + + (unsigned __int128)(v1) + (unsigned __int128)(v2); \ + (hi) = (uint64_t)(fmaz >> 64); \ + (lo) = (uint64_t)fmaz; \ + } while (0) + +/* + * Compute x1*y1+x2*y2+v1+v2. Operands are 64-bit, and result is 128-bit, + * with high word in "hi" and low word in "lo". + * + * Callers should ensure that the two inner products, and the v1 and v2 + * operands, are multiple of 4 (this is not used by this specific definition + * but may help other implementations). + */ +#define FMA2(hi, lo, x1, y1, x2, y2, v1, v2) do { \ + unsigned __int128 fmaz; \ + fmaz = (unsigned __int128)(x1) * (unsigned __int128)(y1) \ + + (unsigned __int128)(x2) * (unsigned __int128)(y2) \ + + (unsigned __int128)(v1) + (unsigned __int128)(v2); \ + (hi) = (uint64_t)(fmaz >> 64); \ + (lo) = (uint64_t)fmaz; \ + } while (0) + +#elif BR_UMUL128 + +#include <intrin.h> + +#define FMA1(hi, lo, x, y, v1, v2) do { \ + uint64_t fmahi, fmalo; \ + unsigned char fmacc; \ + fmalo = _umul128((x), (y), &fmahi); \ + fmacc = _addcarry_u64(0, fmalo, (v1), &fmalo); \ + _addcarry_u64(fmacc, fmahi, 0, &fmahi); \ + fmacc = _addcarry_u64(0, fmalo, (v2), &(lo)); \ + _addcarry_u64(fmacc, fmahi, 0, &(hi)); \ + } while (0) + +/* + * Normally we should use _addcarry_u64() for FMA2 too, but it makes + * Visual Studio crash. Instead we use this version, which leverages + * the fact that the vx operands, and the products, are multiple of 4. + * This is unfortunately slower. + */ +#define FMA2(hi, lo, x1, y1, x2, y2, v1, v2) do { \ + uint64_t fma1hi, fma1lo; \ + uint64_t fma2hi, fma2lo; \ + uint64_t fmatt; \ + fma1lo = _umul128((x1), (y1), &fma1hi); \ + fma2lo = _umul128((x2), (y2), &fma2hi); \ + fmatt = (fma1lo >> 2) + (fma2lo >> 2) \ + + ((v1) >> 2) + ((v2) >> 2); \ + (lo) = fmatt << 2; \ + (hi) = fma1hi + fma2hi + (fmatt >> 62); \ + } while (0) + +/* + * The FMA2 macro definition we would prefer to use, but it triggers + * an internal compiler error in Visual Studio 2015. + * +#define FMA2(hi, lo, x1, y1, x2, y2, v1, v2) do { \ + uint64_t fma1hi, fma1lo; \ + uint64_t fma2hi, fma2lo; \ + unsigned char fmacc; \ + fma1lo = _umul128((x1), (y1), &fma1hi); \ + fma2lo = _umul128((x2), (y2), &fma2hi); \ + fmacc = _addcarry_u64(0, fma1lo, (v1), &fma1lo); \ + _addcarry_u64(fmacc, fma1hi, 0, &fma1hi); \ + fmacc = _addcarry_u64(0, fma2lo, (v2), &fma2lo); \ + _addcarry_u64(fmacc, fma2hi, 0, &fma2hi); \ + fmacc = _addcarry_u64(0, fma1lo, fma2lo, &(lo)); \ + _addcarry_u64(fmacc, fma1hi, fma2hi, &(hi)); \ + } while (0) + */ + +#endif + +#define MASK62 ((uint64_t)0x3FFFFFFFFFFFFFFF) +#define MUL62_lo(x, y) (((uint64_t)(x) * (uint64_t)(y)) & MASK62) + +/* + * Subtract b from a, and return the final carry. If 'ctl32' is 0, then + * a[] is kept unmodified, but the final carry is still computed and + * returned. + */ +static uint32_t +i62_sub(uint64_t *a, const uint64_t *b, size_t num, uint32_t ctl32) +{ + uint64_t cc, mask; + size_t u; + + cc = 0; + ctl32 = -ctl32; + mask = (uint64_t)ctl32 | ((uint64_t)ctl32 << 32); + for (u = 0; u < num; u ++) { + uint64_t aw, bw, dw; + + aw = a[u]; + bw = b[u]; + dw = aw - bw - cc; + cc = dw >> 63; + dw &= MASK62; + a[u] = aw ^ (mask & (dw ^ aw)); + } + return (uint32_t)cc; +} + +/* + * Montgomery multiplication, over arrays of 62-bit values. The + * destination array (d) must be distinct from the other operands + * (x, y and m). All arrays are in little-endian format (least + * significant word comes first) over 'num' words. + */ +static void +montymul(uint64_t *d, const uint64_t *x, const uint64_t *y, + const uint64_t *m, size_t num, uint64_t m0i) +{ + uint64_t dh; + size_t u, num4; + + num4 = 1 + ((num - 1) & ~(size_t)3); + memset(d, 0, num * sizeof *d); + dh = 0; + for (u = 0; u < num; u ++) { + size_t v; + uint64_t f, xu; + uint64_t r, zh; + uint64_t hi, lo; + + xu = x[u] << 2; + f = MUL62_lo(d[0] + MUL62_lo(x[u], y[0]), m0i) << 2; + + FMA2(hi, lo, xu, y[0], f, m[0], d[0] << 2, 0); + r = hi; + + for (v = 1; v < num4; v += 4) { + FMA2(hi, lo, xu, y[v + 0], + f, m[v + 0], d[v + 0] << 2, r << 2); + r = hi + (r >> 62); + d[v - 1] = lo >> 2; + FMA2(hi, lo, xu, y[v + 1], + f, m[v + 1], d[v + 1] << 2, r << 2); + r = hi + (r >> 62); + d[v + 0] = lo >> 2; + FMA2(hi, lo, xu, y[v + 2], + f, m[v + 2], d[v + 2] << 2, r << 2); + r = hi + (r >> 62); + d[v + 1] = lo >> 2; + FMA2(hi, lo, xu, y[v + 3], + f, m[v + 3], d[v + 3] << 2, r << 2); + r = hi + (r >> 62); + d[v + 2] = lo >> 2; + } + for (; v < num; v ++) { + FMA2(hi, lo, xu, y[v], f, m[v], d[v] << 2, r << 2); + r = hi + (r >> 62); + d[v - 1] = lo >> 2; + } + + zh = dh + r; + d[num - 1] = zh & MASK62; + dh = zh >> 62; + } + i62_sub(d, m, num, (uint32_t)dh | NOT(i62_sub(d, m, num, 0))); +} + +/* + * Conversion back from Montgomery representation. + */ +static void +frommonty(uint64_t *x, const uint64_t *m, size_t num, uint64_t m0i) +{ + size_t u, v; + + for (u = 0; u < num; u ++) { + uint64_t f, cc; + + f = MUL62_lo(x[0], m0i) << 2; + cc = 0; + for (v = 0; v < num; v ++) { + uint64_t hi, lo; + + FMA1(hi, lo, f, m[v], x[v] << 2, cc); + cc = hi << 2; + if (v != 0) { + x[v - 1] = lo >> 2; + } + } + x[num - 1] = cc >> 2; + } + i62_sub(x, m, num, NOT(i62_sub(x, m, num, 0))); +} + +/* see inner.h */ +uint32_t +br_i62_modpow_opt(uint32_t *x31, const unsigned char *e, size_t elen, + const uint32_t *m31, uint32_t m0i31, uint64_t *tmp, size_t twlen) +{ + size_t u, mw31num, mw62num; + uint64_t *x, *m, *t1, *t2; + uint64_t m0i; + uint32_t acc; + int win_len, acc_len; + + /* + * Get modulus size, in words. + */ + mw31num = (m31[0] + 31) >> 5; + mw62num = (mw31num + 1) >> 1; + + /* + * In order to apply this function, we must have enough room to + * copy the operand and modulus into the temporary array, along + * with at least two temporaries. If there is not enough room, + * switch to br_i31_modpow(). We also use br_i31_modpow() if the + * modulus length is not at least four words (94 bits or more). + */ + if (mw31num < 4 || (mw62num << 2) > twlen) { + /* + * We assume here that we can split an aligned uint64_t + * into two properly aligned uint32_t. Since both types + * are supposed to have an exact width with no padding, + * then this property must hold. + */ + size_t txlen; + + txlen = mw31num + 1; + if (twlen < txlen) { + return 0; + } + br_i31_modpow(x31, e, elen, m31, m0i31, + (uint32_t *)tmp, (uint32_t *)tmp + txlen); + return 1; + } + + /* + * Convert x to Montgomery representation: this means that + * we replace x with x*2^z mod m, where z is the smallest multiple + * of the word size such that 2^z >= m. We want to reuse the 31-bit + * functions here (for constant-time operation), but we need z + * for a 62-bit word size. + */ + for (u = 0; u < mw62num; u ++) { + br_i31_muladd_small(x31, 0, m31); + br_i31_muladd_small(x31, 0, m31); + } + + /* + * Assemble operands into arrays of 62-bit words. Note that + * all the arrays of 62-bit words that we will handle here + * are without any leading size word. + * + * We also adjust tmp and twlen to account for the words used + * for these extra arrays. + */ + m = tmp; + x = tmp + mw62num; + tmp += (mw62num << 1); + twlen -= (mw62num << 1); + for (u = 0; u < mw31num; u += 2) { + size_t v; + + v = u >> 1; + if ((u + 1) == mw31num) { + m[v] = (uint64_t)m31[u + 1]; + x[v] = (uint64_t)x31[u + 1]; + } else { + m[v] = (uint64_t)m31[u + 1] + + ((uint64_t)m31[u + 2] << 31); + x[v] = (uint64_t)x31[u + 1] + + ((uint64_t)x31[u + 2] << 31); + } + } + + /* + * Compute window size. We support windows up to 5 bits; for a + * window of size k bits, we need 2^k+1 temporaries (for k = 1, + * we use special code that uses only 2 temporaries). + */ + for (win_len = 5; win_len > 1; win_len --) { + if ((((uint32_t)1 << win_len) + 1) * mw62num <= twlen) { + break; + } + } + + t1 = tmp; + t2 = tmp + mw62num; + + /* + * Compute m0i, which is equal to -(1/m0) mod 2^62. We were + * provided with m0i31, which already fulfills this property + * modulo 2^31; the single expression below is then sufficient. + */ + m0i = (uint64_t)m0i31; + m0i = MUL62_lo(m0i, (uint64_t)2 + MUL62_lo(m0i, m[0])); + + /* + * Compute window contents. If the window has size one bit only, + * then t2 is set to x; otherwise, t2[0] is left untouched, and + * t2[k] is set to x^k (for k >= 1). + */ + if (win_len == 1) { + memcpy(t2, x, mw62num * sizeof *x); + } else { + uint64_t *base; + + memcpy(t2 + mw62num, x, mw62num * sizeof *x); + base = t2 + mw62num; + for (u = 2; u < ((unsigned)1 << win_len); u ++) { + montymul(base + mw62num, base, x, m, mw62num, m0i); + base += mw62num; + } + } + + /* + * Set x to 1, in Montgomery representation. We again use the + * 31-bit code. + */ + br_i31_zero(x31, m31[0]); + x31[(m31[0] + 31) >> 5] = 1; + br_i31_muladd_small(x31, 0, m31); + if (mw31num & 1) { + br_i31_muladd_small(x31, 0, m31); + } + for (u = 0; u < mw31num; u += 2) { + size_t v; + + v = u >> 1; + if ((u + 1) == mw31num) { + x[v] = (uint64_t)x31[u + 1]; + } else { + x[v] = (uint64_t)x31[u + 1] + + ((uint64_t)x31[u + 2] << 31); + } + } + + /* + * We process bits from most to least significant. At each + * loop iteration, we have acc_len bits in acc. + */ + acc = 0; + acc_len = 0; + while (acc_len > 0 || elen > 0) { + int i, k; + uint32_t bits; + uint64_t mask1, mask2; + + /* + * Get the next bits. + */ + k = win_len; + if (acc_len < win_len) { + if (elen > 0) { + acc = (acc << 8) | *e ++; + elen --; + acc_len += 8; + } else { + k = acc_len; + } + } + bits = (acc >> (acc_len - k)) & (((uint32_t)1 << k) - 1); + acc_len -= k; + + /* + * We could get exactly k bits. Compute k squarings. + */ + for (i = 0; i < k; i ++) { + montymul(t1, x, x, m, mw62num, m0i); + memcpy(x, t1, mw62num * sizeof *x); + } + + /* + * Window lookup: we want to set t2 to the window + * lookup value, assuming the bits are non-zero. If + * the window length is 1 bit only, then t2 is + * already set; otherwise, we do a constant-time lookup. + */ + if (win_len > 1) { + uint64_t *base; + + memset(t2, 0, mw62num * sizeof *t2); + base = t2 + mw62num; + for (u = 1; u < ((uint32_t)1 << k); u ++) { + uint64_t mask; + size_t v; + + mask = -(uint64_t)EQ(u, bits); + for (v = 0; v < mw62num; v ++) { + t2[v] |= mask & base[v]; + } + base += mw62num; + } + } + + /* + * Multiply with the looked-up value. We keep the product + * only if the exponent bits are not all-zero. + */ + montymul(t1, x, t2, m, mw62num, m0i); + mask1 = -(uint64_t)EQ(bits, 0); + mask2 = ~mask1; + for (u = 0; u < mw62num; u ++) { + x[u] = (mask1 & x[u]) | (mask2 & t1[u]); + } + } + + /* + * Convert back from Montgomery representation. + */ + frommonty(x, m, mw62num, m0i); + + /* + * Convert result into 31-bit words. + */ + for (u = 0; u < mw31num; u += 2) { + uint64_t zw; + + zw = x[u >> 1]; + x31[u + 1] = (uint32_t)zw & 0x7FFFFFFF; + if ((u + 1) < mw31num) { + x31[u + 2] = (uint32_t)(zw >> 31); + } + } + return 1; +} + +#else + +/* see inner.h */ +uint32_t +br_i62_modpow_opt(uint32_t *x31, const unsigned char *e, size_t elen, + const uint32_t *m31, uint32_t m0i31, uint64_t *tmp, size_t twlen) +{ + size_t mwlen; + + mwlen = (m31[0] + 63) >> 5; + if (twlen < mwlen) { + return 0; + } + return br_i31_modpow_opt(x31, e, elen, m31, m0i31, + (uint32_t *)tmp, twlen << 1); +} + +#endif + +/* see inner.h */ +uint32_t +br_i62_modpow_opt_as_i31(uint32_t *x31, const unsigned char *e, size_t elen, + const uint32_t *m31, uint32_t m0i31, uint32_t *tmp, size_t twlen) +{ + /* + * As documented, this function expects the 'tmp' argument to be + * 64-bit aligned. This is OK since this function is internal (it + * is not part of BearSSL's public API). + */ + return br_i62_modpow_opt(x31, e, elen, m31, m0i31, + (uint64_t *)tmp, twlen >> 1); +} |