From 86a4479929bd5f3fa5ce2cabfe6a4ebb53944df4 Mon Sep 17 00:00:00 2001 From: Peter Dettman Date: Fri, 4 Sep 2020 23:57:27 +0700 Subject: 'safegcd' modular inversion --- crypto/src/math/raw/Mod.cs | 554 ++++++++++++++++++++++++++++++++++++++------- crypto/src/math/raw/Nat.cs | 66 ++++++ 2 files changed, 537 insertions(+), 83 deletions(-) (limited to 'crypto/src/math/raw') diff --git a/crypto/src/math/raw/Mod.cs b/crypto/src/math/raw/Mod.cs index 197b5c82b..1a1524c11 100644 --- a/crypto/src/math/raw/Mod.cs +++ b/crypto/src/math/raw/Mod.cs @@ -7,74 +7,190 @@ using Org.BouncyCastle.Utilities; namespace Org.BouncyCastle.Math.Raw { + /* + * Modular inversion as implemented in this class is based on the paper "Fast constant-time gcd + * computation and modular inversion" by Daniel J. Bernstein and Bo-Yin Yang. + */ + internal abstract class Mod { private static readonly SecureRandom RandomSource = new SecureRandom(); - public static void Invert(uint[] p, uint[] x, uint[] z) + private const int M30 = 0x3FFFFFFF; + private const ulong M32UL = 0xFFFFFFFFUL; + + [Obsolete("Will be removed")] + public static void Add(uint[] p, uint[] x, uint[] y, uint[] z) { int len = p.Length; - if (Nat.IsZero(len, x)) - throw new ArgumentException("cannot be 0", "x"); - if (Nat.IsOne(len, x)) + uint c = Nat.Add(len, x, y, z); + if (c != 0) { - Array.Copy(x, 0, z, 0, len); - return; + Nat.SubFrom(len, p, z); } + } - uint[] u = Nat.Copy(len, x); - uint[] a = Nat.Create(len); - a[0] = 1; - int ac = 0; + public static void CheckedModOddInverse(uint[] m, uint[] x, uint[] z) + { + if (0 == ModOddInverse(m, x, z)) + throw new ArithmeticException("Inverse does not exist."); + } - if ((u[0] & 1) == 0) - { - InversionStep(p, u, len, a, ref ac); - } - if (Nat.IsOne(len, u)) + public static void CheckedModOddInverseVar(uint[] m, uint[] x, uint[] z) + { + if (!ModOddInverseVar(m, x, z)) + throw new ArithmeticException("Inverse does not exist."); + } + + public static uint Inverse32(uint d) + { + Debug.Assert((d & 1) == 1); + + //int x = d + (((d + 1) & 4) << 1); // d.x == 1 mod 2**4 + uint x = d; // d.x == 1 mod 2**3 + x *= 2 - d * x; // d.x == 1 mod 2**6 + x *= 2 - d * x; // d.x == 1 mod 2**12 + x *= 2 - d * x; // d.x == 1 mod 2**24 + x *= 2 - d * x; // d.x == 1 mod 2**48 + Debug.Assert(d * x == 1); + return x; + } + + [Obsolete("Use 'CheckedModOddInverseVar' instead")] + public static void Invert(uint[] m, uint[] x, uint[] z) + { + CheckedModOddInverseVar(m, x, z); + } + + public static uint ModOddInverse(uint[] m, uint[] x, uint[] z) + { + int len32 = m.Length; + Debug.Assert(len32 > 0); + Debug.Assert((m[0] & 1) != 0); + Debug.Assert(m[len32 - 1] != 0); + + int bits = (len32 << 5) - Integers.NumberOfLeadingZeros((int)m[len32 - 1]); + int len30 = (bits + 29) / 30; + int m0Inv30x4 = -(int)Inverse32(m[0]) << 2; + + int[] t = new int[4]; + int[] D = new int[len30]; + int[] E = new int[len30]; + int[] F = new int[len30]; + int[] G = new int[len30]; + int[] M = new int[len30]; + + E[0] = 1; + Encode30(bits, x, 0, G, 0); + Encode30(bits, m, 0, M, 0); + Array.Copy(M, 0, F, 0, len30); + + int eta = -1; + int maxDivsteps = GetMaximumDivsteps(bits); + + for (int divSteps = 0; divSteps < maxDivsteps; divSteps += 30) { - InversionResult(p, ac, a, z); - return; + eta = Divsteps30(eta, F[0], G[0], t); + UpdateDE30(len30, D, E, t, m0Inv30x4, M); + UpdateFG30(len30, F, G, t); } - uint[] v = Nat.Copy(len, p); - uint[] b = Nat.Create(len); - int bc = 0; + int signF = F[len30 - 1] >> 31; + Debug.Assert(-1 == signF | 0 == signF); - int uvLen = len; + CNegate30(len30, signF, F); + CNegate30(len30, signF, D); - for (;;) + Decode30(bits, D, 0, z, 0); + + int signD = D[len30 - 1] >> 31; + Debug.Assert(-1 == signD | 0 == signD); + + signD += (int)Nat.CAdd(len32, signD, z, m, z); + Debug.Assert(0 == signD & 0 != Nat.LessThan(len32, z, m)); + + return (uint)(EqualTo(len30, F, 1) & EqualToZero(len30, G)); + } + + public static bool ModOddInverseVar(uint[] m, uint[] x, uint[] z) + { + int len32 = m.Length; + Debug.Assert(len32 > 0); + Debug.Assert((m[0] & 1) != 0); + Debug.Assert(m[len32 - 1] != 0); + + int bits = (len32 << 5) - Integers.NumberOfLeadingZeros((int)m[len32 - 1]); + int len30 = (bits + 29) / 30; + int m0Inv30x4 = -(int)Inverse32(m[0]) << 2; + + int[] t = new int[4]; + int[] D = new int[len30]; + int[] E = new int[len30]; + int[] F = new int[len30]; + int[] G = new int[len30]; + int[] M = new int[len30]; + + E[0] = 1; + Encode30(bits, x, 0, G, 0); + Encode30(bits, m, 0, M, 0); + Array.Copy(M, 0, F, 0, len30); + + int clzG = Integers.NumberOfLeadingZeros(G[len30 - 1] | 1) - (len30 * 30 + 2 - bits); + int eta = -1 - clzG; + int lenDE = len30, lenFG = len30; + int maxDivsteps = GetMaximumDivsteps(bits); + + int divsteps = 0; + while (!IsZero(lenFG, G)) { - while (u[uvLen - 1] == 0 && v[uvLen - 1] == 0) - { - --uvLen; - } + if (divsteps >= maxDivsteps) + return false; - if (Nat.Gte(len, u, v)) - { - Nat.SubFrom(len, v, u); - Debug.Assert((u[0] & 1) == 0); - ac += Nat.SubFrom(len, b, a) - bc; - InversionStep(p, u, uvLen, a, ref ac); - if (Nat.IsOne(len, u)) - { - InversionResult(p, ac, a, z); - return; - } - } - else + divsteps += 30; + + eta = Divsteps30Var(eta, F[0], G[0], t); + UpdateDE30(lenDE, D, E, t, m0Inv30x4, M); + UpdateFG30(lenFG, F, G, t); + + int fn = F[lenFG - 1]; + int gn = G[lenFG - 1]; + + int cond = (lenFG - 2) >> 31; + cond |= fn ^ (fn >> 31); + cond |= gn ^ (gn >> 31); + + if (cond == 0) { - Nat.SubFrom(len, u, v); - Debug.Assert((v[0] & 1) == 0); - bc += Nat.SubFrom(len, a, b) - ac; - InversionStep(p, v, uvLen, b, ref bc); - if (Nat.IsOne(len, v)) - { - InversionResult(p, bc, b, z); - return; - } + F[lenFG - 2] |= fn << 30; + G[lenFG - 2] |= gn << 30; + --lenFG; } } + + int signF = F[lenFG - 1] >> 31; + Debug.Assert(-1 == signF | 0 == signF); + + if (0 != signF) + { + Negate30(lenFG, F); + Negate30(lenDE, D); + } + + if (!IsOne(lenFG, F)) + return false; + + Decode30(bits, D, 0, z, 0); + + int signD = D[lenDE - 1] >> 31; + Debug.Assert(-1 == signD | 0 == signD); + + if (signD < 0) + { + signD += (int)Nat.AddTo(len32, m, z); + } + Debug.Assert(0 == signD && !Nat.Gte(len32, z, m)); + + return true; } public static uint[] Random(uint[] p) @@ -101,74 +217,346 @@ namespace Org.BouncyCastle.Math.Raw return s; } - public static void Add(uint[] p, uint[] x, uint[] y, uint[] z) + [Obsolete("Will be removed")] + public static void Subtract(uint[] p, uint[] x, uint[] y, uint[] z) { int len = p.Length; - uint c = Nat.Add(len, x, y, z); + int c = Nat.Sub(len, x, y, z); if (c != 0) { - Nat.SubFrom(len, p, z); + Nat.AddTo(len, p, z); } } - public static void Subtract(uint[] p, uint[] x, uint[] y, uint[] z) + private static void CNegate30(int len, int cond, int[] D) { - int len = p.Length; - int c = Nat.Sub(len, x, y, z); - if (c != 0) + Debug.Assert(len > 0); + Debug.Assert(D.Length >= len); + + int last = len - 1; + long cd = 0L; + + for (int i = 0; i < last; ++i) { - Nat.AddTo(len, p, z); + cd += (D[i] ^ cond) - cond; + D[i] = (int)cd & M30; cd >>= 30; } + + cd += (D[last] ^ cond) - cond; + D[last] = (int)cd; } - private static void InversionResult(uint[] p, int ac, uint[] a, uint[] z) + private static void Decode30(int bits, int[] x, int xOff, uint[] z, int zOff) { - if (ac < 0) + Debug.Assert(bits > 0); + + int avail = 0; + ulong data = 0L; + + while (bits > 0) { - Nat.Add(p.Length, a, p, z); + while (avail < System.Math.Min(32, bits)) + { + data |= (ulong)x[xOff++] << avail; + avail += 30; + } + + z[zOff++] = (uint)data; data >>= 32; + avail -= 32; + bits -= 32; } - else + } + + private static int Divsteps30(int eta, int f0, int g0, int[] t) + { + int u = 1, v = 0, q = 0, r = 1; + int f = f0, g = g0; + + for (int i = 0; i < 30; ++i) { - Array.Copy(a, 0, z, 0, p.Length); + Debug.Assert((f & 1) == 1); + Debug.Assert((u * f0 + v * g0) == f << i); + Debug.Assert((q * f0 + r * g0) == g << i); + + int p = -(g & 1); + int s = eta >> 31; + + int c1 = p & s; + int c2 = p & ~s; + + eta = (eta ^ c1) - (c1 + 1); + + g += f & c2; + q += u & c2; + r += v & c2; + + g -= f & c1; + q -= u & c1; + r -= v & c1; + + f += g & c1; + u += q & c1; + v += r & c1; + + g >>= 1; + u <<= 1; + v <<= 1; } + + t[0] = u; + t[1] = v; + t[2] = q; + t[3] = r; + + return eta; } - private static void InversionStep(uint[] p, uint[] u, int uLen, uint[] x, ref int xc) + private static int Divsteps30Var(int eta, int f0, int g0, int[] t) { - int len = p.Length; - int count = 0; - while (u[0] == 0) + int u = 1, v = 0, q = 0, r = 1; + int f = f0, g = g0, m, w, x, y, z; + int i = 30, limit, zeros; + + for (;;) { - Nat.ShiftDownWord(uLen, u, 0); - count += 32; + // Use a sentinel bit to count zeros only up to i. + zeros = Integers.NumberOfTrailingZeros(g | (-1 << i)); + + g >>= zeros; + u <<= zeros; + v <<= zeros; + eta -= zeros; + i -= zeros; + + if (i <= 0) + break; + + Debug.Assert((f & 1) == 1); + Debug.Assert((g & 1) == 1); + Debug.Assert((u * f0 + v * g0) == f << (30 - i)); + Debug.Assert((q * f0 + r * g0) == g << (30 - i)); + + if (eta < 0) + { + eta = -eta; + x = f; f = g; g = -x; + y = u; u = q; q = -y; + z = v; v = r; r = -z; + + // Handle up to 6 divsteps at once, subject to eta and i. + limit = (eta + 1) > i ? i : (eta + 1); + m = (int)((uint.MaxValue >> (32 - limit)) & 63U); + + w = (f * g * (f * f - 2)) & m; + } + else + { + // Handle up to 4 divsteps at once, subject to eta and i. + limit = (eta + 1) > i ? i : (eta + 1); + m = (int)((uint.MaxValue >> (32 - limit)) & 15U); + + w = f + (((f + 1) & 4) << 1); + w = (-w * g) & m; + } + + g += f * w; + q += u * w; + r += v * w; + + Debug.Assert((g & m) == 0); } + t[0] = u; + t[1] = v; + t[2] = q; + t[3] = r; + + return eta; + } + + private static void Encode30(int bits, uint[] x, int xOff, int[] z, int zOff) + { + Debug.Assert(bits > 0); + + int avail = 0; + ulong data = 0UL; + + while (bits > 0) { - int zeroes = Integers.NumberOfTrailingZeros((int)u[0]); - if (zeroes > 0) + if (avail < System.Math.Min(30, bits)) { - Nat.ShiftDownBits(uLen, u, zeroes, 0); - count += zeroes; + data |= (x[xOff++] & M32UL) << avail; + avail += 32; } + + z[zOff++] = (int)data & M30; data >>= 30; + avail -= 30; + bits -= 30; + } + } + + private static int EqualTo(int len, int[] x, int y) + { + int d = x[0] ^ y; + for (int i = 1; i < len; ++i) + { + d |= x[i]; + } + d = (int)((uint)d >> 1) | (d & 1); + return (d - 1) >> 31; + } + + private static int EqualToZero(int len, int[] x) + { + int d = 0; + for (int i = 0; i < len; ++i) + { + d |= x[i]; } + d = (int)((uint)d >> 1) | (d & 1); + return (d - 1) >> 31; + } + + private static int GetMaximumDivsteps(int bits) + { + return (49 * bits + (bits < 46 ? 80 : 47)) / 17; + } - for (int i = 0; i < count; ++i) + private static bool IsOne(int len, int[] x) + { + if (x[0] != 1) { - if ((x[0] & 1) != 0) + return false; + } + for (int i = 1; i < len; ++i) + { + if (x[i] != 0) { - if (xc < 0) - { - xc += (int)Nat.AddTo(len, p, x); - } - else - { - xc += Nat.SubFrom(len, p, x); - } + return false; } + } + return true; + } - Debug.Assert(xc == 0 || xc == -1); - Nat.ShiftDownBit(len, x, (uint)xc); + private static bool IsZero(int len, int[] x) + { + if (x[0] != 0) + { + return false; + } + for (int i = 1; i < len; ++i) + { + if (x[i] != 0) + { + return false; + } + } + return true; + } + + private static void Negate30(int len, int[] D) + { + Debug.Assert(len > 0); + Debug.Assert(D.Length >= len); + + int last = len - 1; + long cd = 0L; + + for (int i = 0; i < last; ++i) + { + cd -= D[i]; + D[i] = (int)cd & M30; cd >>= 30; + } + + cd -= D[last]; + D[last] = (int)cd; + } + + private static void UpdateDE30(int len, int[] D, int[] E, int[] t, int m0Inv30x4, int[] M) + { + Debug.Assert(len > 0); + Debug.Assert(D.Length >= len); + Debug.Assert(E.Length >= len); + Debug.Assert(M.Length >= len); + Debug.Assert(m0Inv30x4 * M[0] == -1 << 2); + + int u = t[0], v = t[1], q = t[2], r = t[3]; + int di, ei, i, md, me; + long cd, ce; + + di = D[0]; + ei = E[0]; + + cd = (long)u * di + (long)v * ei; + ce = (long)q * di + (long)r * ei; + + md = (m0Inv30x4 * (int)cd) >> 2; + me = (m0Inv30x4 * (int)ce) >> 2; + + cd += (long)M[0] * md; + ce += (long)M[0] * me; + + Debug.Assert(((int)cd & M30) == 0); + Debug.Assert(((int)ce & M30) == 0); + + cd >>= 30; + ce >>= 30; + + for (i = 1; i < len; ++i) + { + di = D[i]; + ei = E[i]; + + cd += (long)u * di + (long)v * ei; + ce += (long)q * di + (long)r * ei; + + cd += (long)M[i] * md; + ce += (long)M[i] * me; + + D[i - 1] = (int)cd & M30; cd >>= 30; + E[i - 1] = (int)ce & M30; ce >>= 30; + } + + D[len - 1] = (int)cd; + E[len - 1] = (int)ce; + } + + private static void UpdateFG30(int len, int[] F, int[] G, int[] t) + { + Debug.Assert(len > 0); + Debug.Assert(F.Length >= len); + Debug.Assert(G.Length >= len); + + int u = t[0], v = t[1], q = t[2], r = t[3]; + int fi, gi, i; + long cf, cg; + + fi = F[0]; + gi = G[0]; + + cf = (long)u * fi + (long)v * gi; + cg = (long)q * fi + (long)r * gi; + + Debug.Assert(((int)cf & M30) == 0); + Debug.Assert(((int)cg & M30) == 0); + + cf >>= 30; + cg >>= 30; + + for (i = 1; i < len; ++i) + { + fi = F[i]; + gi = G[i]; + + cf += (long)u * fi + (long)v * gi; + cg += (long)q * fi + (long)r * gi; + + F[i - 1] = (int)cf & M30; cf >>= 30; + G[i - 1] = (int)cg & M30; cg >>= 30; } + + F[len - 1] = (int)cf; + G[len - 1] = (int)cg; } } } diff --git a/crypto/src/math/raw/Nat.cs b/crypto/src/math/raw/Nat.cs index 69942661f..9f2412580 100644 --- a/crypto/src/math/raw/Nat.cs +++ b/crypto/src/math/raw/Nat.cs @@ -447,6 +447,72 @@ namespace Org.BouncyCastle.Math.Raw return true; } + public static uint EqualTo(int len, uint[] x, uint y) + { + uint d = x[0] ^ y; + for (int i = 1; i < len; ++i) + { + d |= x[i]; + } + d = (d >> 1) | (d & 1); + return (uint)(((int)d - 1) >> 31); + } + + public static uint EqualTo(int len, uint[] x, int xOff, uint y) + { + uint d = x[xOff] ^ y; + for (int i = 1; i < len; ++i) + { + d |= x[xOff + i]; + } + d = (d >> 1) | (d & 1); + return (uint)(((int)d - 1) >> 31); + } + + public static uint EqualTo(int len, uint[] x, uint[] y) + { + uint d = 0; + for (int i = 0; i < len; ++i) + { + d |= x[i] ^ y[i]; + } + d = (d >> 1) | (d & 1); + return (uint)(((int)d - 1) >> 31); + } + + public static uint EqualTo(int len, uint[] x, int xOff, uint[] y, int yOff) + { + uint d = 0; + for (int i = 0; i < len; ++i) + { + d |= x[xOff + i] ^ y[yOff + i]; + } + d = (d >> 1) | (d & 1); + return (uint)(((int)d - 1) >> 31); + } + + public static uint EqualToZero(int len, uint[] x) + { + uint d = 0; + for (int i = 0; i < len; ++i) + { + d |= x[i]; + } + d = (d >> 1) | (d & 1); + return (uint)(((int)d - 1) >> 31); + } + + public static uint EqualToZero(int len, uint[] x, int xOff) + { + uint d = 0; + for (int i = 0; i < len; ++i) + { + d |= x[xOff + i]; + } + d = (d >> 1) | (d & 1); + return (uint)(((int)d - 1) >> 31); + } + public static uint[] FromBigInteger(int bits, BigInteger x) { if (x.SignValue < 0 || x.BitLength > bits) -- cgit 1.4.1