| ... | ... | @@ -3,6 +3,8 @@ |
| 3 | 3 | //! Unlike `std.math.big`, these integers have a fixed maximum length and are only designed to be used for modular arithmetic. |
| 4 | 4 | //! Arithmetic operations are meant to run in constant-time for a given modulus, making them suitable for cryptography. |
| 5 | 5 | //! |
| 6 | //! Functions with `Public` in their name are the exception. |
| 7 | //! |
| 6 | 8 | //! Parts of that code was ported from the BSD-licensed crypto/internal/bigmod/nat.go file in the Go language, itself inspired from BearSSL. |
| 7 | 9 | |
| 8 | 10 | const std = @import("std"); |
| ... | ... | @@ -51,8 +53,10 @@ pub const FieldElementError = error{NonCanonical}; |
| 51 | 53 | /// Invalid representation (Montgomery vs non-Montgomery domain.) |
| 52 | 54 | pub const RepresentationError = error{UnexpectedRepresentation}; |
| 53 | 55 | |
| 56 | pub const DivisionByZeroError = error{DivisionByZero}; |
| 57 | |
| 54 | 58 | /// The set of all possible errors `std.crypto.ff` functions can return. |
| 55 | | pub const Error = OverflowError || InvalidModulusError || NullExponentError || FieldElementError || RepresentationError; |
| 59 | pub const Error = OverflowError || InvalidModulusError || NullExponentError || FieldElementError || RepresentationError || DivisionByZeroError; |
| 56 | 60 | |
| 57 | 61 | /// An unsigned big integer with a fixed maximum size (`max_bits`), suitable for cryptographic operations. |
| 58 | 62 | /// Storage rounds up to whole limbs and can hold up to `capacity_bits` bits. |
| ... | ... | @@ -240,6 +244,110 @@ pub fn Uint(comptime max_bits: comptime_int) type { |
| 240 | 244 | return x.conditionalSubWithOverflow(true, y); |
| 241 | 245 | } |
| 242 | 246 | |
| 247 | /// Adds the product of `x` and `y` to `acc`, wrapping at the active width. |
| 248 | /// Returns 1 if any part of the result was discarded. |
| 249 | pub fn mulAddWithOverflow(acc: *Self, x: Self, y: Self) u1 { |
| 250 | assert(x.limbs_len == acc.limbs_len); |
| 251 | assert(y.limbs_len == acc.limbs_len); |
| 252 | const n = acc.limbs_len; |
| 253 | var wide: [2 * max_limbs_count]Limb = undefined; |
| 254 | @memcpy(wide[0..n], acc.limbsConst()); |
| 255 | mulAddWide(wide[0 .. 2 * n], x.limbsConst(), y.limbsConst()); |
| 256 | @memcpy(acc.limbs(), wide[0..n]); |
| 257 | return @intFromBool(!ct.eql(orLimbs(wide[n..][0..n]), 0)); |
| 258 | } |
| 259 | |
| 260 | /// Returns the bit length, or 0 for zero. |
| 261 | pub fn bitLenPublic(x: Self) usize { |
| 262 | var i = x.limbs_len; |
| 263 | while (i != 0) { |
| 264 | i -= 1; |
| 265 | const limb = x.limbsConst()[i]; |
| 266 | if (limb != 0) { |
| 267 | return i * t_bits + t_bits - @clz(@as(TLimb, @intCast(limb))); |
| 268 | } |
| 269 | } |
| 270 | return 0; |
| 271 | } |
| 272 | |
| 273 | /// Returns the number of trailing zero bits, or the active width for zero. |
| 274 | pub fn trailingZeroBitsPublic(x: Self) usize { |
| 275 | for (x.limbsConst(), 0..) |limb, i| { |
| 276 | if (limb != 0) { |
| 277 | return i * t_bits + @ctz(@as(TLimb, @intCast(limb))); |
| 278 | } |
| 279 | } |
| 280 | return x.limbs_len * t_bits; |
| 281 | } |
| 282 | |
| 283 | pub fn shiftRightPublic(x: *Self, shift: usize) void { |
| 284 | const limb_shift = shift / t_bits; |
| 285 | const bit_shift = shift % t_bits; |
| 286 | const x_limbs = x.limbs(); |
| 287 | if (limb_shift >= x_limbs.len) { |
| 288 | @memset(x_limbs, 0); |
| 289 | return; |
| 290 | } |
| 291 | const active = x_limbs.len - limb_shift; |
| 292 | for (0..active) |i| { |
| 293 | const lo = math.shr(Limb, x_limbs[i + limb_shift], bit_shift); |
| 294 | const hi = if (i + limb_shift + 1 < x_limbs.len) |
| 295 | math.shl(Limb, x_limbs[i + limb_shift + 1], t_bits - bit_shift) |
| 296 | else |
| 297 | 0; |
| 298 | x_limbs[i] = @as(TLimb, @truncate(lo | hi)); |
| 299 | } |
| 300 | @memset(x_limbs[active..], 0); |
| 301 | } |
| 302 | |
| 303 | /// Divides by `divisor` in place and returns the remainder. |
| 304 | /// Returns `error.DivisionByZero` without changing the integer if `divisor` is zero. |
| 305 | pub fn divRemPublic(x: *Self, divisor: usize) DivisionByZeroError!usize { |
| 306 | if (divisor == 0) return error.DivisionByZero; |
| 307 | const Wide = @Int(.unsigned, 2 * @bitSizeOf(Limb)); |
| 308 | const x_limbs = x.limbs(); |
| 309 | var rem: Limb = 0; |
| 310 | var i = x.limbs_len; |
| 311 | while (i != 0) { |
| 312 | i -= 1; |
| 313 | const num = (@as(Wide, rem) << t_bits) | x_limbs[i]; |
| 314 | x_limbs[i] = @intCast(num / divisor); |
| 315 | rem = @intCast(num % divisor); |
| 316 | } |
| 317 | return rem; |
| 318 | } |
| 319 | |
| 320 | /// Returns the greatest common divisor, with gcd(x, 0) = x. |
| 321 | /// Returns `error.DivisionByZero` if both operands are zero. |
| 322 | pub fn gcdPublic(x: Self, y: Self) DivisionByZeroError!Self { |
| 323 | if (x.isZero() and y.isZero()) return error.DivisionByZero; |
| 324 | if (x.isZero()) return y; |
| 325 | if (y.isZero()) return x; |
| 326 | |
| 327 | var a = x; |
| 328 | var b = y; |
| 329 | const len = @max(a.limbs_len, b.limbs_len); |
| 330 | a.expandTo(len); |
| 331 | b.expandTo(len); |
| 332 | |
| 333 | const a_shift = a.trailingZeroBitsPublic(); |
| 334 | const b_shift = b.trailingZeroBitsPublic(); |
| 335 | a.shiftRightPublic(a_shift); |
| 336 | b.shiftRightPublic(b_shift); |
| 337 | while (true) { |
| 338 | switch (a.compare(b)) { |
| 339 | .eq => break, |
| 340 | .lt => mem.swap(Self, &a, &b), |
| 341 | .gt => {}, |
| 342 | } |
| 343 | _ = a.subWithOverflow(b); |
| 344 | a.shiftRightPublic(a.trailingZeroBitsPublic()); |
| 345 | } |
| 346 | |
| 347 | a.shiftLeft(@min(a_shift, b_shift)); // GCD always fits in either input, overflow is never an issue |
| 348 | return a; |
| 349 | } |
| 350 | |
| 243 | 351 | fn expandTo(x: *Self, new_len: usize) void { |
| 244 | 352 | assert(new_len >= x.limbs_len and new_len <= x.limbs_buffer.len); |
| 245 | 353 | @memset(x.limbs_buffer[x.limbs_len..new_len], 0); |
| ... | ... | @@ -276,6 +384,40 @@ pub fn Uint(comptime max_bits: comptime_int) type { |
| 276 | 384 | } |
| 277 | 385 | return borrow; |
| 278 | 386 | } |
| 387 | |
| 388 | fn shiftLeft(x: *Self, shift: usize) void { |
| 389 | assert(x.bitLenPublic() + shift <= x.limbs_len * t_bits); |
| 390 | const limb_shift = shift / t_bits; |
| 391 | const bit_shift = shift % t_bits; |
| 392 | const x_limbs = x.limbs(); |
| 393 | var i = x_limbs.len; |
| 394 | while (i != 0) { |
| 395 | i -= 1; |
| 396 | const hi = if (i >= limb_shift) |
| 397 | math.shl(Limb, x_limbs[i - limb_shift], bit_shift) |
| 398 | else |
| 399 | 0; |
| 400 | const lo = if (i >= limb_shift + 1) |
| 401 | math.shr(Limb, x_limbs[i - limb_shift - 1], t_bits - bit_shift) |
| 402 | else |
| 403 | 0; |
| 404 | x_limbs[i] = @as(TLimb, @truncate(hi | lo)); |
| 405 | } |
| 406 | } |
| 407 | |
| 408 | // Shifts in `carry` at the top and returns the low bit shifted out. |
| 409 | fn shiftRightByOneWithCarry(x: *Self, carry: u1) u1 { |
| 410 | var c: Limb = carry; |
| 411 | var i = x.limbs_len; |
| 412 | const x_limbs = x.limbs(); |
| 413 | while (i != 0) { |
| 414 | i -= 1; |
| 415 | const limb = x_limbs[i]; |
| 416 | x_limbs[i] = (limb >> 1) | (c << (t_bits - 1)); |
| 417 | c = @as(u1, @truncate(limb)); |
| 418 | } |
| 419 | return @truncate(c); |
| 420 | } |
| 279 | 421 | }; |
| 280 | 422 | } |
| 281 | 423 | |
| ... | ... | @@ -402,6 +544,30 @@ fn orLimbs(limbs: []const Limb) Limb { |
| 402 | 544 | return t; |
| 403 | 545 | } |
| 404 | 546 | |
| 547 | // Adds `x * y` to `z` and returns the carry. |
| 548 | fn mulAddLimb(z: []Limb, x: []const Limb, y: Limb) Limb { |
| 549 | assert(z.len == x.len); |
| 550 | var carry: Limb = 0; |
| 551 | for (z, x) |*z_limb, x_limb| { |
| 552 | const wide = ct.mulWide(x_limb, y); |
| 553 | var z_lo = @addWithOverflow(z_limb.*, wide.lo); |
| 554 | var z_hi = wide.hi +% z_lo[1]; |
| 555 | z_lo = @addWithOverflow(z_lo[0], carry); |
| 556 | z_hi +%= z_lo[1]; |
| 557 | z_limb.* = @as(TLimb, @truncate(z_lo[0])); |
| 558 | carry = (z_hi << 1) | (z_lo[0] >> t_bits); |
| 559 | } |
| 560 | return carry; |
| 561 | } |
| 562 | |
| 563 | // Adds `x * y` to the low `x.len` limbs of `z`. The upper limbs need no initialization. |
| 564 | fn mulAddWide(z: []Limb, x: []const Limb, y: []const Limb) void { |
| 565 | assert(z.len == x.len + y.len); |
| 566 | for (y, 0..) |y_limb, i| { |
| 567 | z[i + x.len] = mulAddLimb(z[i..][0..x.len], x, y_limb); |
| 568 | } |
| 569 | } |
| 570 | |
| 405 | 571 | /// A modulus, defining a finite field. |
| 406 | 572 | /// All operations within the field are performed modulo this modulus, without heap allocations. |
| 407 | 573 | /// `max_bits` represents the number of bits in the maximum value the modulus can be set to. |
| ... | ... | @@ -1219,6 +1385,70 @@ test "field element decoding" { |
| 1219 | 1385 | } |
| 1220 | 1386 | } |
| 1221 | 1387 | |
| 1388 | test "Uint bit measurement and shifts" { |
| 1389 | if (builtin.zig_backend == .stage2_c) return error.SkipZigTest; |
| 1390 | |
| 1391 | const U = Uint(256); |
| 1392 | try testing.expectEqual(0, U.zero.bitLenPublic()); |
| 1393 | try testing.expectEqual(U.max_limbs_count * t_bits, U.zero.trailingZeroBitsPublic()); |
| 1394 | try testing.expectEqual(t_bits, U.zero.normalize().trailingZeroBitsPublic()); |
| 1395 | for ([_]usize{ 0, t_bits - 1, t_bits, t_bits + 1, 255 }) |offset| { |
| 1396 | const x = try U.fromPrimitive(u256, math.shl(u256, 1, offset)); |
| 1397 | try testing.expectEqual(offset + 1, x.bitLenPublic()); |
| 1398 | try testing.expectEqual(offset, x.trailingZeroBitsPublic()); |
| 1399 | } |
| 1400 | |
| 1401 | const v: u256 = (1 << 255) | (1 << t_bits) | 3; |
| 1402 | for ([_]usize{ 0, 1, t_bits - 1, t_bits, t_bits + 1, U.max_limbs_count * t_bits, 10_000 }) |shift| { |
| 1403 | var x = try U.fromPrimitive(u256, v); |
| 1404 | x.shiftRightPublic(shift); |
| 1405 | try testing.expectEqual(math.shr(u256, v, shift), try x.toPrimitive(u256)); |
| 1406 | try testing.expectEqual(U.max_limbs_count, x.limbs_len); |
| 1407 | try expectWellFormedLimbs(x); |
| 1408 | } |
| 1409 | for ([_]usize{ 0, 1, t_bits - 1, t_bits, t_bits + 1 }) |shift| { |
| 1410 | var x = try U.fromPrimitive(u256, (1 << t_bits) + 3); |
| 1411 | x.shiftLeft(shift); |
| 1412 | try testing.expectEqual(math.shl(u256, (1 << t_bits) + 3, shift), try x.toPrimitive(u256)); |
| 1413 | } |
| 1414 | |
| 1415 | var x = (try U.fromPrimitive(u8, 5)).normalize(); |
| 1416 | try testing.expectEqual(1, x.shiftRightByOneWithCarry(1)); |
| 1417 | try testing.expectEqual((1 << (t_bits - 1)) + 2, try x.toPrimitive(u64)); |
| 1418 | try testing.expectEqual(0, x.shiftRightByOneWithCarry(0)); |
| 1419 | x.shiftRightPublic(t_bits); |
| 1420 | try testing.expect(x.isZero()); |
| 1421 | try testing.expectEqual(1, x.limbs_len); |
| 1422 | } |
| 1423 | |
| 1424 | test "Uint division and gcd" { |
| 1425 | if (builtin.zig_backend == .stage2_c) return error.SkipZigTest; |
| 1426 | |
| 1427 | const U = Uint(256); |
| 1428 | const v: u256 = (1 << 255) + (1 << t_bits) + 13; |
| 1429 | for ([_]usize{ 1, 7, math.maxInt(TLimb), (1 << t_bits) + 5, math.maxInt(usize) }) |divisor| { |
| 1430 | var x = try U.fromPrimitive(u256, v); |
| 1431 | try testing.expectEqual(v % divisor, try x.divRemPublic(divisor)); |
| 1432 | try testing.expectEqual(v / divisor, try x.toPrimitive(u256)); |
| 1433 | try expectWellFormedLimbs(x); |
| 1434 | } |
| 1435 | var x = (try U.fromPrimitive(u16, 1000)).normalize(); |
| 1436 | try testing.expectError(error.DivisionByZero, x.divRemPublic(0)); |
| 1437 | try testing.expectEqual(6, try x.divRemPublic(7)); |
| 1438 | try testing.expectEqual(142, try x.toPrimitive(u16)); |
| 1439 | try testing.expectEqual(1, x.limbs_len); |
| 1440 | try testing.expectEqual(142, try x.divRemPublic(1000)); |
| 1441 | try testing.expectEqual(0, try x.divRemPublic(7)); |
| 1442 | |
| 1443 | for ([_][3]u256{ .{ 0, 5, 5 }, .{ 6, 6, 6 }, .{ 17, 4, 1 }, .{ 15 << 200, 21 << 100, 3 << 100 } }) |c| { |
| 1444 | const a = (try U.fromPrimitive(u256, c[0])).normalize(); |
| 1445 | const b = try U.fromPrimitive(u256, c[1]); |
| 1446 | try testing.expectEqual(c[2], try (try a.gcdPublic(b)).toPrimitive(u256)); |
| 1447 | try testing.expectEqual(c[2], try (try b.gcdPublic(a)).toPrimitive(u256)); |
| 1448 | } |
| 1449 | try testing.expectError(error.DivisionByZero, U.zero.gcdPublic(U.zero)); |
| 1450 | } |
| 1451 | |
| 1222 | 1452 | test "Uint addition and multiply-add" { |
| 1223 | 1453 | if (builtin.zig_backend == .stage2_c) return error.SkipZigTest; |
| 1224 | 1454 | |
| ... | ... | @@ -1236,4 +1466,29 @@ test "Uint addition and multiply-add" { |
| 1236 | 1466 | try testing.expect(wrapped.eql(full)); |
| 1237 | 1467 | try testing.expectEqual(1, wrapped.addWithOverflow(one)); |
| 1238 | 1468 | try testing.expect(wrapped.isZero()); |
| 1469 | |
| 1470 | for ([_]struct { acc: U, x: U, y: U, overflow: u1, expected: U }{ |
| 1471 | .{ .acc = U.zero, .x = full, .y = one, .overflow = 0, .expected = full }, |
| 1472 | .{ .acc = full, .x = U.zero, .y = full, .overflow = 0, .expected = full }, |
| 1473 | .{ .acc = full, .x = full, .y = U.zero, .overflow = 0, .expected = full }, |
| 1474 | .{ .acc = U.zero, .x = full, .y = full, .overflow = 1, .expected = one }, |
| 1475 | .{ .acc = full, .x = one, .y = one, .overflow = 1, .expected = U.zero }, |
| 1476 | }) |c| { |
| 1477 | var acc = c.acc; |
| 1478 | try testing.expectEqual(c.overflow, acc.mulAddWithOverflow(c.x, c.y)); |
| 1479 | try testing.expect(acc.eql(c.expected)); |
| 1480 | try expectWellFormedLimbs(acc); |
| 1481 | } |
| 1482 | const a: u256 = (1 << 100) + 3; |
| 1483 | const b: u256 = (1 << 80) + 7; |
| 1484 | var acc = try U.fromPrimitive(u256, a); |
| 1485 | try testing.expectEqual(0, acc.mulAddWithOverflow(try U.fromPrimitive(u256, a), try U.fromPrimitive(u256, b))); |
| 1486 | try testing.expectEqual(a + a * b, try acc.toPrimitive(u256)); |
| 1487 | |
| 1488 | var short = one.normalize(); |
| 1489 | const high = (try U.fromPrimitive(u64, 1 << (t_bits - 1))).normalize(); |
| 1490 | const two = (try U.fromPrimitive(u8, 2)).normalize(); |
| 1491 | try testing.expectEqual(1, short.mulAddWithOverflow(high, two)); |
| 1492 | try testing.expectEqual(1, short.limbs_len); |
| 1493 | try testing.expect(short.isOne()); |
| 1239 | 1494 | } |