|
| 1 | +// This Source Code Form is subject to the terms of the Mozilla Public |
| 2 | +// License, v. 2.0. If a copy of the MPL was not distributed with this |
| 3 | +// file, You can obtain one at https://mozilla.org/MPL/2.0/. |
| 4 | + |
| 5 | +//! Implementation of Finite Field GF(2^8) aka GF(256) |
| 6 | +//! |
| 7 | +//! We use the Rijndael (AES) polynomial (x^8 + x^4 + x^3 + x + 1) as our |
| 8 | +//! irreducible polynomial. |
| 9 | +//! |
| 10 | +//! We only implement enough operations so that we can implement |
| 11 | +//! shamir secret sharing. This keeps the surface area small and |
| 12 | +//! minimizes the amount of auditing we need to do. |
| 13 | +//! |
| 14 | +//! For a basic overview of galois fields, the docs from the |
| 15 | +//! [gf256 crate](https://docs.rs/gf256/0.3.0/gf256/gf/index.html) |
| 16 | +//! are excellent. |
| 17 | +//! |
| 18 | +//! For a more comprehensive introduction to the math, these |
| 19 | +//! [lecture notes](https://web.stanford.edu/~marykw/classes/CS250_W19/readings/Forney_Introduction_to_Finite_Fields.pdf) |
| 20 | +//! are handy. |
| 21 | +
|
| 22 | +// Don't tell me what operations to use in my implementations |
| 23 | +#![expect(clippy::suspicious_arithmetic_impl)] |
| 24 | + |
| 25 | +use core::fmt::{self, Binary, Display, Formatter, LowerHex, UpperHex}; |
| 26 | +use core::ops::{Add, AddAssign, Div, Mul, MulAssign, Sub}; |
| 27 | +use rand::Rng; |
| 28 | +use rand::distributions::{Distribution, Standard}; |
| 29 | +use subtle::ConstantTimeEq; |
| 30 | +use zeroize::Zeroize; |
| 31 | + |
| 32 | +/// An element in a finite field of prime power 2^8 |
| 33 | +/// |
| 34 | +/// We explicitly don't enable the equality operators to prevent ourselves from |
| 35 | +/// accidentally using those instead of the constant time ones. |
| 36 | +#[repr(transparent)] |
| 37 | +#[derive(Debug, Clone, Copy, Zeroize)] |
| 38 | +pub struct Gf256(u8); |
| 39 | + |
| 40 | +impl Gf256 { |
| 41 | + pub fn new(n: u8) -> Gf256 { |
| 42 | + Gf256(n) |
| 43 | + } |
| 44 | + |
| 45 | + /// Return the underlying u8. |
| 46 | + pub fn into_u8(self) -> u8 { |
| 47 | + self.0 |
| 48 | + } |
| 49 | + |
| 50 | + /// Return the multiplicative inverse (`self^-1`) of self |
| 51 | + /// |
| 52 | + /// ```text |
| 53 | + /// self * self^-1 = 1 |
| 54 | + /// ``` |
| 55 | + /// |
| 56 | + /// By Fermat's little theorem: `self^-1 = self^254 for GF(2^8)`. |
| 57 | + /// We calculate `self^254` in a simple, unrolled fashion. |
| 58 | + /// |
| 59 | + /// This strategy was borrowed from <https://github.com/dsprenkels/sss/blob/16c3fdb175497b25eb90b966991fa7ff19fbdcfe/hazmat.c#L247-L266> |
| 60 | + #[rustfmt::skip] |
| 61 | + pub fn invert(&self) -> Gf256 { |
| 62 | + let mut result = *self * *self; // self^2 |
| 63 | + let temp_4 = result * result; // self^4 |
| 64 | + result = temp_4 * temp_4; // self^8 |
| 65 | + let temp_9 = result * *self; // self^9 |
| 66 | + result *= result; // self^16; |
| 67 | + result *= temp_9; // self^25; |
| 68 | + let temp_50 = result * result; // self^50; |
| 69 | + result = temp_50 * temp_50; // self^100; |
| 70 | + result *= result; // self^200; |
| 71 | + result = result * temp_50 * temp_4; // self^254 |
| 72 | + |
| 73 | + result |
| 74 | + } |
| 75 | +} |
| 76 | + |
| 77 | +pub const ZERO: Gf256 = Gf256(0); |
| 78 | +pub const ONE: Gf256 = Gf256(1); |
| 79 | + |
| 80 | +impl Display for Gf256 { |
| 81 | + fn fmt(&self, f: &mut Formatter<'_>) -> std::fmt::Result { |
| 82 | + write!(f, "{}", self.0) |
| 83 | + } |
| 84 | +} |
| 85 | + |
| 86 | +impl LowerHex for Gf256 { |
| 87 | + fn fmt(&self, f: &mut Formatter<'_>) -> fmt::Result { |
| 88 | + write!(f, "{:02x}", self.0) |
| 89 | + } |
| 90 | +} |
| 91 | + |
| 92 | +impl UpperHex for Gf256 { |
| 93 | + fn fmt(&self, f: &mut Formatter<'_>) -> fmt::Result { |
| 94 | + write!(f, "{:02X}", self.0) |
| 95 | + } |
| 96 | +} |
| 97 | + |
| 98 | +impl Binary for Gf256 { |
| 99 | + fn fmt(&self, f: &mut Formatter<'_>) -> fmt::Result { |
| 100 | + write!(f, "{:08b}", self.0) |
| 101 | + } |
| 102 | +} |
| 103 | + |
| 104 | +impl Distribution<Gf256> for Standard { |
| 105 | + fn sample<R: Rng + ?Sized>(&self, rng: &mut R) -> Gf256 { |
| 106 | + Gf256(rng.r#gen()) |
| 107 | + } |
| 108 | +} |
| 109 | + |
| 110 | +impl ConstantTimeEq for Gf256 { |
| 111 | + fn ct_eq(&self, other: &Self) -> subtle::Choice { |
| 112 | + self.0.ct_eq(&other.0) |
| 113 | + } |
| 114 | +} |
| 115 | + |
| 116 | +impl Add for Gf256 { |
| 117 | + type Output = Self; |
| 118 | + |
| 119 | + /// This is a fundamental operation. |
| 120 | + /// |
| 121 | + /// In `GF(2^n)`, we always do carryless addition `mod 2` which ends up |
| 122 | + /// being exactly equal to bitwise xor. |
| 123 | + fn add(self, rhs: Self) -> Self::Output { |
| 124 | + Gf256(self.0 ^ rhs.0) |
| 125 | + } |
| 126 | +} |
| 127 | + |
| 128 | +impl Add<&Gf256> for Gf256 { |
| 129 | + type Output = Self; |
| 130 | + fn add(self, rhs: &Gf256) -> Self::Output { |
| 131 | + self + *rhs |
| 132 | + } |
| 133 | +} |
| 134 | + |
| 135 | +impl AddAssign for Gf256 { |
| 136 | + fn add_assign(&mut self, rhs: Self) { |
| 137 | + *self = *self + rhs |
| 138 | + } |
| 139 | +} |
| 140 | + |
| 141 | +impl MulAssign for Gf256 { |
| 142 | + fn mul_assign(&mut self, rhs: Self) { |
| 143 | + *self = *self * rhs |
| 144 | + } |
| 145 | +} |
| 146 | + |
| 147 | +impl Sub for Gf256 { |
| 148 | + type Output = Self; |
| 149 | + |
| 150 | + /// Subtraction is identical to addition in `GF(2^n)` |
| 151 | + fn sub(self, rhs: Self) -> Self { |
| 152 | + Gf256(self.0 ^ rhs.0) |
| 153 | + } |
| 154 | +} |
| 155 | + |
| 156 | +impl Div for Gf256 { |
| 157 | + type Output = Self; |
| 158 | + |
| 159 | + fn div(self, rhs: Self) -> Self::Output { |
| 160 | + self * rhs.invert() |
| 161 | + } |
| 162 | +} |
| 163 | + |
| 164 | +/// This is an efficient carryless multiplication + reduction algorithm. |
| 165 | +/// |
| 166 | +/// It is essentially long multiplication of polynomial elements in GF(2^8) |
| 167 | +/// represented as unsigned bytes (u8). For every bit set in `rhs` we multiply |
| 168 | +/// our original value by `x^i` where i is the bit position of `rhs` and `x` is |
| 169 | +/// our placeholder variable for GF(2^8) operations, and then add the result to |
| 170 | +/// an accumulator initialized to zero. If our accumulator exceeds 255, meaning |
| 171 | +/// it is no longer inside the finite field, we reduce it modulo our irreducible |
| 172 | +/// (prime) polynomial: |
| 173 | +/// |
| 174 | +/// ```text |
| 175 | +/// m(x) = x^8 + x^4 + x^3 + x + 1 |
| 176 | +/// ``` |
| 177 | +/// |
| 178 | +/// As an example, with our accumulator named `product`: |
| 179 | +/// |
| 180 | +/// ```text |
| 181 | +/// product = 0 |
| 182 | +/// self = 0x21 = 0b0010_0001 = x^5 + 1 |
| 183 | +/// rhs = 0x12 = 0b0001_0010 = x^4 + x |
| 184 | +/// |
| 185 | +/// step 1 = ((x^5 + 1) * x) mod m(x) |
| 186 | +/// = x^6 + x |
| 187 | +/// |
| 188 | +/// product += step1 |
| 189 | +/// |
| 190 | +/// step 2 = ((x^5 + 1) * x^4) mod m(x) |
| 191 | +/// = (x^9 + x^4) mod (x^8 + x^4 + x^3 +x + 1) |
| 192 | +/// = x^5 + x^2 + x |
| 193 | +/// |
| 194 | +/// product += step2 |
| 195 | +/// = x^6 + x + x^5 + x^2 + x |
| 196 | +/// = x^6 + x^5 + x^2 |
| 197 | +/// = 0b0110_0100 |
| 198 | +/// = 0x64 |
| 199 | +/// ``` |
| 200 | +/// |
| 201 | +/// See the `test_docs_example` unit test at the bottom of this file to confirm |
| 202 | +/// this math. |
| 203 | +/// |
| 204 | +/// In this algorithm, we use the value `0x1b`, which is our irreducible |
| 205 | +/// polynomial without the high term. We do this because of the following |
| 206 | +/// equality in GF(2^8), which follows from long-division of polynomials: |
| 207 | +/// |
| 208 | +/// ```text |
| 209 | +/// `x^8 mod m(x) = x^4 + x^3 + x + 1 = 0x1b` |
| 210 | +/// ``` |
| 211 | +/// |
| 212 | +/// A rationale and description of this algorithm can be found in sections 7.9 |
| 213 | +/// and 7.10 of the lecture notes at |
| 214 | +/// <https://engineering.purdue.edu/kak/compsec/NewLectures/Lecture7.pdf> |
| 215 | +/// |
| 216 | +/// It is also explained on a wikipedia section about the |
| 217 | +/// [AES finite field](https://en.wikipedia.org/wiki/Finite_field_arithmetic#Rijndael's_(AES)_finite_field) |
| 218 | +/// |
| 219 | +/// This is roughly the same algorithm used in |
| 220 | +/// [vsss-rs](https://github.com/mikelodder7/vsss-rs) and many other existing |
| 221 | +/// implementations. |
| 222 | +impl Mul for Gf256 { |
| 223 | + type Output = Self; |
| 224 | + fn mul(mut self, mut rhs: Self) -> Self::Output { |
| 225 | + let mut product = 0u8; |
| 226 | + for _ in 0..8 { |
| 227 | + // If the LSB of rhs is 1, then add `self` to `product`. |
| 228 | + product ^= 0x0u8.wrapping_sub(rhs.0 & 1) & self.0; |
| 229 | + |
| 230 | + // Divide rhs polynomial by `x`. |
| 231 | + // Note that we just utilized this bit we right shifted away in the |
| 232 | + // previous line. |
| 233 | + rhs.0 >>= 1; |
| 234 | + |
| 235 | + // Track if the high bit is currently set in `self`. If it is, then |
| 236 | + // we have an `x^7` term, and multiplying by `x` will require a |
| 237 | + // modulo reduction to stay within our field. |
| 238 | + let carry: u8 = self.0 >> 7; |
| 239 | + |
| 240 | + // Shift `self` left a bit. |
| 241 | + // This multiplies `self` by `x`. |
| 242 | + self.0 <<= 1; |
| 243 | + |
| 244 | + // If there was a carry, then we need to add `x^8 mod m(x)`, which |
| 245 | + // equals `0x1b`, to `self`. |
| 246 | + self.0 ^= 0x0u8.wrapping_sub(carry) & 0x1b; |
| 247 | + |
| 248 | + // Note that we haven't actually added `self` to `product` yet, as |
| 249 | + // we don't know if the corresponding bit in `rhs` was set. |
| 250 | + // That happens the next time through the loop. |
| 251 | + } |
| 252 | + |
| 253 | + Gf256(product) |
| 254 | + } |
| 255 | +} |
| 256 | + |
| 257 | +impl Mul<&Gf256> for Gf256 { |
| 258 | + type Output = Gf256; |
| 259 | + fn mul(self, rhs: &Self) -> Self::Output { |
| 260 | + self * *rhs |
| 261 | + } |
| 262 | +} |
| 263 | + |
| 264 | +#[cfg(test)] |
| 265 | +mod tests { |
| 266 | + use super::*; |
| 267 | + |
| 268 | + #[test] |
| 269 | + fn test_all_multiplication() { |
| 270 | + // This file contains the output of all products of (0..255)x(0..255) |
| 271 | + // as produced by a snapshot of vsss_rs version 5.10 and validated via |
| 272 | + // gf256 version 0.3.0. |
| 273 | + let products_tabbed: &str = include_str!(concat!( |
| 274 | + env!("CARGO_MANIFEST_DIR"), |
| 275 | + "/test_input/gf256_multiplication_test_vector.txt" |
| 276 | + )); |
| 277 | + |
| 278 | + let mut products = products_tabbed.split("\t"); |
| 279 | + for i in 0..=255 { |
| 280 | + for j in 0..=255 { |
| 281 | + let a = Gf256(i); |
| 282 | + let b = Gf256(j); |
| 283 | + let c = a * b; |
| 284 | + assert_eq!( |
| 285 | + c.0, |
| 286 | + products.next().unwrap().parse::<u8>().unwrap() |
| 287 | + ); |
| 288 | + } |
| 289 | + } |
| 290 | + } |
| 291 | + |
| 292 | + #[test] |
| 293 | + fn test_all_multiplicative_inverses() { |
| 294 | + // Can't divide by zero |
| 295 | + for i in 1..=255u8 { |
| 296 | + let a = Gf256::new(i); |
| 297 | + assert_eq!((a * a.invert()).ct_eq(&ONE).unwrap_u8(), 1); |
| 298 | + |
| 299 | + // Division is the same as multiplying by the inverse |
| 300 | + assert_eq!((a / a).ct_eq(&ONE).unwrap_u8(), 1); |
| 301 | + } |
| 302 | + } |
| 303 | + |
| 304 | + #[test] |
| 305 | + fn test_docs_example() { |
| 306 | + let a = Gf256::new(0x21); |
| 307 | + let b = Gf256::new(0x12); |
| 308 | + |
| 309 | + let expected = Gf256::new(0x64); |
| 310 | + let product = a * b; |
| 311 | + |
| 312 | + println!("product = {:#x?}", product.0); |
| 313 | + |
| 314 | + assert_eq!(expected.ct_eq(&product).unwrap_u8(), 1); |
| 315 | + } |
| 316 | +} |
0 commit comments