Coverage for src/sbk/gf_util.py : 99%
Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
1# This file is part of the sbk project
2# https://github.com/mbarkhau/sbk
3#
4# Copyright (c) 2019-2021 Manuel Barkhau (mbarkhau@gmail.com) - MIT License
5# SPDX-License-Identifier: MIT
7"""Galois Field arithmetic functions."""
9import typing as typ
11from . import gf_lut
13# https://en.wikipedia.org/wiki/Finite_field_arithmetic#Rijndael's_finite_field
14#
15# x**8 + x**4 + x**3 + x + 1 (0b100011011 = 0x11B)
16#
17# Rijndael uses the characteristic 2 finite field with 256 elements, which can
18# also be called the Galois field GF(2**8). It employs the following reducing
19# polynomial for multiplication:
20#
21# For example, 0x53 * 0xCA = 0x01 in Rijndael's field because
22#
23# 0x53 = 0b01010011
24# 0xCA = 0b11001010
25#
26# (x6 + x4 + x + 1)(x7 + x6 + x3 + x)
27# = (x13 + x12 + x9 + x7) + (x11 + x10 + x7 + x5) + (x8 + x7 + x4 + x2) + (x7 + x6 + x3 + x)
28# = x13 + x12 + x9 + x11 + x10 + x5 + x8 + x4 + x2 + x6 + x3 + x
29# = x13 + x12 + x11 + x10 + x9 + x8 + x6 + x5 + x4 + x3 + x2 + x
30# and
31#
32# x13 + x12 + x11 + x10 + x9 + x8 + x6 + x5 + x4 + x3 + x2 + x
33# mod x8 + x4 + x3 + x1 + 1
34# = (0b11111101111110 mod 100011011) = 0x3F7E mod 0x011B = 0x0001
35#
36# , which can be demonstrated through long division (shown using binary
37# notation, since it lends itself well to the task. Notice that exclusive OR
38# is applied in the example and not arithmetic subtraction, as one might use
39# in grade-school long division.):
40#
41# 11111101111110 (mod) 100011011
42# ^10001101100000
43# 1110000011110
44# ^1000110110000
45# 110110101110
46# ^100011011000
47# 10101110110
48# ^10001101100
49# 100011010
50# ^100011011
51# 1
53# The multiplicative inverse for an element a of a finite field can be calculated
54# a number of different ways:
55#
56# Since the nonzero elements of GF(p^n) form a finite group with respect to multiplication,
57# a^((p^n)−1) = 1 (for a != 0), thus the inverse of a is a^((p^n)−2).
60RIJNDAEL_REDUCING_POLYNOMIAL = 0x011B
63def div_slow(a: int, b: int) -> int:
64 # long division
65 val = a
66 divisor = b
67 assert divisor > 0
69 while divisor < val:
70 divisor = divisor << 1
72 mask = 1
73 while mask < divisor:
74 mask = mask << 1
75 mask = mask >> 1
77 while divisor > 0xFF:
78 if (val & mask) > 0:
79 val = val ^ divisor
81 divisor = divisor >> 1
82 mask = mask >> 1
84 return val
87DIV_LUT: typ.Dict[int, int] = {}
90def div(a: int, b: int) -> int:
91 assert 0 <= a < 65536, a
92 assert b > 0, b
94 key = a * 65536 + b
95 if key not in DIV_LUT:
96 DIV_LUT[key] = div_slow(a, b)
98 return DIV_LUT[key]
101def div_by_rrp(val: int) -> int:
102 return div(val, RIJNDAEL_REDUCING_POLYNOMIAL)
105def mul_slow(a: int, b: int) -> int:
106 res = 0
107 while a > 0:
108 if a & 1 != 0:
109 res = res ^ b
110 a = a // 2
111 b = b * 2
113 return div_by_rrp(res)
116def mul(a: int, b: int) -> int:
117 assert 0 <= a < 256, a
118 assert 0 <= b < 256, b
120 mul_lut = gf_lut.MUL_LUT
121 if not mul_lut:
122 gf_lut.init_mul_lut()
124 return mul_lut[a][b]
127def pow_slow(a: int, b: int) -> int:
128 res = 1
129 n = b
130 while n > 0:
131 res = mul(res, a)
132 n -= 1
133 return res
136def inverse_slow(val: int) -> int:
137 """Calculate multiplicative inverse in GF(256).
139 Since the nonzero elements of GF(p^n) form a finite group with
140 respect to multiplication,
142 a^((p^n)−1) = 1 (for a != 0)
144 thus the inverse of a is
146 a^((p^n)−2).
147 """
148 if val == 0:
149 return 0
151 exp = 2 ** 8 - 2
152 inv = pow_slow(val, exp)
153 assert mul(val, inv) == 1
154 return inv
157def inverse(val: int) -> int:
158 return gf_lut.MUL_INVERSE_LUT[val]
161# The Euclidean GCD algorithm is based on the principle that the
162# greatest common divisor of two numbers does not change if the larger
163# number is replaced by its difference with the smaller number. For
164# example,
165#
166# GCD(252) == 21 # 252 = 21 × 12
167# GCD(105) == 21 # 105 = 21 × 5
168#
169# also
170#
171# GCD(252 − 105) == 21
172# GCD(147) == 21
173#
174# Since this replacement reduces the larger of the two numbers,
175# repeating this process gives successively smaller pairs of numbers
176# until the two numbers become equal. When that occurs, they are the
177# GCD of the original two numbers.
178#
179# By reversing the steps, the GCD can be expressed as a sum of the two
180# original numbers each multiplied by a positive or negative integer,
181# e.g., 21 = 5 × 105 + (−2) × 252. The fact that the GCD can always be
182# expressed in this way is known as Bézout's identity.
185class XGCDResult(typ.NamedTuple):
186 g: int
187 s: int # sometimes called x
188 t: int # sometimes called y
191def xgcd(a: int, b: int) -> XGCDResult:
192 """Extended euclidien greatest common denominator."""
193 # pylint: disable=unpacking-non-sequence
195 if a == 0:
196 return XGCDResult(b, 0, 1)
198 g, s, t = xgcd(b % a, a)
200 q = b // a
201 res = XGCDResult(g=g, s=t - q * s, t=s)
202 assert res.s * a + res.t * b == res.g
203 return res