GF2p8 Кодирование Рида-Соломона в изоморфном представлении поля

Коды Рида-Соломона считаются через экспоненты и логарифмы по таблице. А у нас с внедренеием инструкций Intel GFNI появилось две другие возможности: использовать инструкцию умножения и аффинные преобразования.

Для примера я взял коды РС из проекта библиотеки генеации QR-кодов. Там используется полином поля G(2^8) P = 0x11D. Кодирование производится в цикле. При переходе к операциям умножения в поле 0x11D, цикл кодрования выглядит так.

	for(i = 0; i < data_length; i++) {
		uint8_t fb =  data[i] ^ Z[0];
		for(j = 0; j < ecc_length-1; j++) {
			Z[j] = Z[j+1] ^ gf2p8_mul(fb, g[j], 0x11D);
		}
		Z[j] = gf2p8_mul(fb, g[j], 0x11D);
	}
	for(i = 0; i < ecc_length; i++) {
		ecc[i] = Z[i];
	}

Преобразование алгоритма для расчета в изоморфном представлении поля с образующим полиномом 0x11B:

	const uint64_t M = 0xFFAACC88F0A0C080;
	const uint64_t Mr= 0xFFAACC88F0A0C080;

	for(i = 0; i < data_length; i++) {
		uint8_t fb =  affine(M,data[i]) ^ Z[0];
		for(j = 0; j < ecc_length-1; j++) {
			Z[j] = Z[j+1] ^ gf2p8_mul(fb, g[j], 0x11B);
		}
		Z[j] = gf2p8_mul(fb, g[j], 0x11B);
	}
	for(i = 0; i < ecc_length; i++) {
		ecc[i] = affine(Mr, Z[i]);
	}

Данное преобразование предполагает так же изоморфные преобразования компонент вектора генератора кодов (g). Векторизация алгоритма может быть выполнена по вложенному циклу.

[1]The comeback of Reed Solomon codes N.Drucker, Sh.Gueron, V.Krasnov

Приложение

Довеском идет код для расчета остатка от деления на 255. Когда в контроллере нет полиномиального умножения, приходится использовать способ расчета умножения и инверсии с использованием таблиц логарифмов и экспонент.

// Multiplication in Galua Fields
int gf_mul(int x, int y) {
    if (x == 0 || y == 0) return 0;
    return exp_[(log_[x] + log_[y])%255];
}
// Inversion in Galua Fields
int gf_inv(int x) {
    return exp_[255 - log_[x]];
}
// Division in Galua Fields
int gf_div(int x, int y) {
    if (x == 0) return 0;
    return exp_[(log_[x] + 255 - log_[y]) % 255];
}

Самая неудобная операция в этом случае - вычисление остатка от деления на 255. Избавился от деления.

// Остаток от деления на 255
uint8_t mod255(unsigned v) {
	uint32_t q = (v*0x1010102UL)>>(24);
	return q;
}

// Остаток от деления на 255
uint8_t mod255(unsigned v) {
	return v + (v>>8);
}

А этот вариант самый простой, применим для сложения логарифмов.

Такая вот оптимизация.

[2] Исходники на Github

GF2p8 Нахождение корней образующего полинома

Понадобилось больше теории, чтобы перейти к композитным полям. Но все попадаются сатьи индусов с ведической алгеброй или индусов из исследовательского центра IBM в Индии. В единственном нормальном учебнике, который мне попался по теории Галуа и кодам Рида-Соломона, слово "корень" употребляется всего один раз, вообще один раз на все сто-пятьсот страниц, и без определения.

Допустим A -- один из корней образующего полинома p(x) = x⁸+x⁴+x³+x+1. При сложении коэффициентов при степенях используется арифметика по модулю 2. Степени A образуют полиномиальный базис [1,A,A²,A³,A⁴,A⁵,A⁶,A⁷], который далее может использоваться для вычисления матриц изоморфных преобразований и представления элементов поля.

Что значит корень? Это значит, что подстановка "корня" в образующий поле полином (8-ого порядка) дает ноль:
P(A) = с₈A⁸ + с₇A⁷ + ... + c₂A² + c₁A¹ + c₀=0.
Здесь коэффициенты при степенях (c_k) те же, что и в образующем полиноме.
Как найти корень? Чтобы найти в одно действие какой-то корень надо отбросить от образующего полинома степень 8, получим полином 7й степени – это один из корней (x⁴+x³+x+1). Отбрасывание старшего порядка равносильно операции A = x⁸ mod P(x). Остальные корни находим методом возведения полученного корня в квадрат, рекурсивно. Всего будет 8 корней для полинома 8 порядка.

// Проверка на корень
int gf2p8_is_root(uint8_t a, uint32_t Poly)
{
	uint8_t v=1, r=0;
	uint32_t p = Poly;
	while (p) {
		if(p&1) r ^= v;
		v = gf2p8_mul(v, a, Poly);
		p>>=1;
	}
	return (r==0);
}

Полином поля можно представить иначе, с учетом знания о корнях, как произведение:
P(x)= (x-A)(x-A²)(x-A⁴)(x-A⁸)(x-A¹⁶)(x-A³²)(x-A⁶⁴)(x-A¹²⁸).

Все корни образующего полинома поля GF(2^m)образуют циклическую группу по операции возведения в квадрат.

Для полинома восьмого порядка сразу можно выделить четыре корня: A =x¹, x², x⁴, x⁸ mod P(x). Остальные четыре корня надо считать.

Свойства поля Галуа GF(p^m):

1. Существует порождающий примитивный элемент поля (геренатор поля), такой что любой не нулевой элемент поля может быть выражен, как степень порождающего элемента поля. Мультипликативная группа циклична.
2. Для любого элемента (g) поля GF(p^m) справедливо: g^pm-1 = 1, g^-1 = g^pm-2.
3. Любой элемент поля может быть выражен суммой из m элементов поля. Для выражения элементов поля используются m элементов, составляющих полиномиальный базис. Элемент поля можно представить через полиномиальный базис: a₀ +a₁A +a₂A² +...+a_m−1A^m−1.
4. В поле GF(p^m) выполняется равенство: (a + b)^p = a^p + b^p (mod p)
5. Если элемент A поля GF(p^m) - корень нередуцируемого полинома P(x) степени m, то остальными корнями P(x) будут элементы: A^p, A^p2, . . . , A^pm-1.
6. Для любого простого числа p и любого нередуцируемого полинома P(x) степени m существует только одно поле Галуа GF(p^m). Поля Галуа GF(p^m), образованные различными нередуцируемыми полиномами степени m изоморфны.

Ниже привожу все корни от 30 нередуцируемых полиномов 8-ого порядка, и корни для полиномов 4-ого порядка.

Таблица корней полинома 8-ого порядка
11B: 02 04 10 1B 5E E4 4D FA
11D: 02 04 10 1D 4C 9D 5F 85
12B: 02 04 10 2B E9 EE FB C1
12D: 02 04 10 2D E5 5D F6 75
139: 02 04 10 39 9C 4E 8F 72
13F: 02 04 10 3F 96 91 84 AA
14D: 02 04 10 4D F8 E3 EB AB
15F: 02 04 10 5F 86 9A 95 C0
163: 02 04 10 63 F3 F4 E1 93
165: 02 04 10 65 8B 8C 99 ED
169: 02 04 10 69 03 05 11 68
171: 02 04 10 71 F3 2F E0 5B
177: 02 04 10 77 8B 3A 98 48
17B: 02 04 10 7B 03 05 11 7A
187: 02 04 10 87 6F 1D D6 34
18B: 02 04 10 8B AA 1A CF E3
18D: 02 04 10 8D 3D 7F 60 B8
19F: 02 04 10 9F 7D 59 94 38
1A3: 02 04 10 A3 70 CD 4A 43
1A9: 02 04 10 A9 18 E9 CD 82
1B1: 02 04 10 B1 66 5F B8 27
1BD: 02 04 10 BD 95 E6 A1 78
1C3: 02 04 10 C3 1F 96 A0 FD
1CF: 02 04 10 CF 76 A0 58 56
1D7: 02 04 10 D7 EB 9E CC 79
1DD: 02 04 10 DD 61 4D D3 35
1E7: 02 04 10 E7 7D 4B EA 2C
1F3: 02 04 10 F3 C1 1C A3 9A
1F5: 02 04 10 F5 28 7E 67 D3
1F9: 02 04 10 F9 6F 1B BC 26

Таблица корней полинома 4-ого порядка
13: 2 4 3 5
19: 2 4 9 E
1F: 2 4 F 8

GF2p8 Аффинные преобразования в Композитные поля GF((2^4)^2)

Расширение и сужение поля

Идея в том, чтобы вместо операции над 8-и битными значениям использовать операции с 4-х битными значениями. Как осуществить переход?

В теории полей Галуа есть понятие расширение поля. Расширением GF(2^4) образованного полиномом q(x) до GF(2^8) c образующим полиномом p(z) является полином второго порядка GF((2^4)^2) s(y), но действующий с 4х битными коэффициентами. Запишем так:
s(y) = y²+t·y+n. – полином расширения поля, второго порядка,
где t – “trap” n – “norm” коэффициенты 4 бит. Коэффициенты хотелось бы выбрать таким образом, чтобы один из них (t) обратился в единицу.

Необходимо описать, как выглядят операции в поле: умножение и инверсия. Все числа в надстройке поля представляются парами: старшая и младшая часть:
a = a₁·y+a₀
y - символизирует сдвиг на 4 разряда.

Умножение в композитном поле выглядит следующим образом:
(a₁·y+a₀)(b₁·y+b₀) = (a₁·b₁)·y² + (a₀·b₁+b₀·a₁)·y + a₀·b₀

Используя редуцирование по полиному y2+t·y+n =0, получаем:
(a₁·y+a₀)(b₁·y+b₀) = (a₀·b₁ + b₀·a₁ + a₁·b₁·t)·y + (a₀·b₀ + a₁·b₁·n) = r₁·y + r₀
Результат представлен в виде пары чисел разрядностью 4 бита.

Таким образом, можно выразить умножение, через умножение нижнем поле GF(2^4):
r₁ = (a₀·b₁ + b₀·a₁ + a₁·b₁·t) mod q(x)
r₀ = (a₀·b₀ + a₁·b₁·n) mod q(x)

Также можем выразить операцию инверсии, через инверсию в нижнем поле:
(a₁·y + a₀)^-1 = (d·a₁)·y + d·(a₁·t + a₀),
где d = (a₁²·n + a₁·a₀·t + a₀²)^-1 –- инверсия в поле GF(2^4).
r₁ = (d·a₁) mod q(x)
r₀ = d·(a₁·t + a₀) mod q(x)

Можно проверить справедливость формулы инверсии с использованием операции умножения.

Надо определить коэффициенты {t, n}, выбрать полином в поле GF(2^4).
Всего существует 4 значения n удовлетворяющих требованию t=1 и три нередуцируемх полинома четвертой степени. Для каждого набора параметров можно составить 8 матриц трансформации.

Нередуцируемые полиномы 4-ого порядка:
13: q(x) =x4+x+1
19: q(x) =x4+x3+1
1F: q(x) =x4+x3+x2+x+1

Вариант расчета: t=1, n=1100 – используется в литературе, применительно к полям GF(2^8)/0x11B (AES) и в сочетании с полиномом GF(2^4)/0x13.

Нас интересует составление матриц изоморфного преобразования между полем GF(2^8) и композитным полем GF((2^4)^2)

Для реализации нам понадобятся следующие алгоритмы:
1. умножение и инверсия в поле расширения;
2. умножение и инверсия в поле основания GF(2^4);
3. выбор генератора поля в композитном поле GF((2^4)^2);
4. составление и проверка матрицы трансформации между полем GF(2^8) и композитным полем GF((2^4)^2).

Операции в поле GF(2^4) можно реализовать методом табличной подстановки. Ниже привожу таблицы инверсии и редуцирования. Кроме того, можно реализовать таблицу квадратов и таблицу умножения на константу. Любую функцию одного переменного мы можем записать LUT таблицей (LookUp Table), четыре входа. При описаний логики для FPGA, четырех-битные LUT – оптимальный случай, поскольку, логика на четыре входа и таблица на четыре входа умещается в один логический элемент. При разработке алгоритмов с использованием векторных инструкций выигрыш в том, что вся таблица помещается в вектор 16 байт и для подстановки можно использовать функцию SHUFFLE.

Таблицы инверсии для полиномов 4-ого порядка:
13: 0,1,9,E,D,B,7,6,F,2,C,5,A,4,3,8
19: 0,1,C,8,6,F,4,E,3,D,B,A,2,9,7,5
1F: 0,1,F,A,8,6,5,9,4,7,3,E,D,C,B,2
Таблицы редуцирования старшей части полиномиального умножения:
13: 0,3,6,5,C,F,A,9,B,8,D,E,7,4,1,2
19: 0,9,B,2,F,6,4,D,7,E,C,5,8,1,3,A
1f: 0,F,1,E,2,D,3,C,4,B,5,A,6,9,7,8

Выбор генератора композитного поля можно совместить с генерацией таблицы инверсии, как было показано ранее. В данном случае нам нужно установить соответствие между (a) примитивным элементом поля GF(2^8) и (b) примитивным элементом поля GF((2^4)^2).

Алгоритм подбора генератора поля GF((2^4)^2):
1. Выбираем следующий примитивный элемент (b) поля GF((2^4)^2).
2. Составляем каким-то чудом матрицу прямого и обратного преобразования.
3. Выполняем проверку: произведение прямой и обратной матрицы дает единичную. 
Проверяем действительно ли полученное преобразование – изоморфизм. 
Если условие не выполняется, выбираем следующий примитивный элемент (b),
возвращаемся к шагу 1.

Свойства изоморфизма
    map(a×b) = map(a)×map(b)
    map(a⊕b) = map(a)⊕map(b)

Проверить можно гомоморфизм групповой по отношению к операциям сложения и умножения, для i,j принадлежащих множеству [0..2^k-1] должно выполняться:
a^t = aⁱ+a^j => b^t = bⁱ+b^j
a^i+j = aⁱ·a^j => b^i+j = bⁱ·b^j

Выбор генератора поля позволяет составить таблицу инверсии в композитном поле. Критерием того, что элемент поля может использоваться в качестве генератора, является целостность таблицы, степени генератора поля позволяют получить все элементы поля за исключением нуля.

Multiplicative Inverse GF((2^4)^2) 0x13, N=1100
 00 01 09 0E 0D 0B 07 06 0F 02 0C 05 0A 04 03 08
 AA A0 C7 CB 52 57 4F 4B FE F1 3E 3D 82 8A 2D 2F
 55 E9 50 E7 6A 8F 6C 87 B1 45 BA 41 2E 1E 2C 1F
 66 FA F5 60 96 48 4C 9F 5E D8 D5 5B 3F 1B 1A 3C
 BB 2B 7D C9 B0 29 7A C5 35 97 4E 17 36 9E 4A 16
 22 E6 14 56 E8 20 53 15 A9 75 D4 3B 72 A3 38 D9
 33 F4 E5 B7 EB BC 30 FB D3 9C 24 8E 26 86 DE 95
 99 B4 5C A2 A8 59 BF 90 D1 EF 46 C4 C8 42 E1 DC
 CC DA 1C 8B AF 9A 6D 27 C0 D7 1D 83 A5 93 6B 25
 77 BE A4 8D DF 6F 34 49 B5 70 85 AE 69 D2 4D 37
 11 AB 73 5D 92 8C B3 F2 74 58 10 A1 B8 FD 9B 84
 44 28 F3 A6 71 98 E4 63 AC FC 2A 40 65 EA 91 76
 88 D6 E2 F8 7B 47 CA 12 7C 43 C6 13 80 DB EC F7
 EE 78 9D 68 5A 3A C1 89 39 5F 81 CD 7F E0 6E 94
 DD 7E C2 F9 B6 62 51 23 54 21 BD 64 CE F6 D0 79
 FF 19 A7 B2 61 32 ED CF C3 E3 31 67 B9 AD 18 F0

Алгоритм (2) поиска матрицы прямого и обратного преобразования привожу ниже. Идея сводится к тому, чтобы из степеней генератора поля (a) составить диагональную матрицу и применив те же степени к генератору поля (b) сформировать по столбцам матрицу изоморфного преобразования. Для составления обратной матрицы поля меняются местами.

/* Алгоритм составления матрицы изоморфного преобразования
    a - примитивный элемент поля А
    P1 - образующий полином поля A
    b - примитивный элемент поля B 
    P2 - образующий полином поля B
 */
uint64_t map_ab(uint8_t a, uint32_t P1, uint8_t b, uint32_t P2)
{
	uint64_t m=0;
	int i;
	for(i=0;i<8;i++){
		int k=1;
		uint8_t v=1, r=1;
		while (v != (1<<(i))) {
			v = gf2p8_mul(v, a, P1);
			r = gf2p8_mul(r, b, P2);
		}
		m |= (uint64_t)r<<(8*(7-i));
	}
	return m_transpose(m);
}

В представленном алгоритме при выполнении операции умножения для композитного поля применяется арифметика композитного поля GF((2^4)^2), для поля AES арифметика поля GF(2^8), функции gf2p8_mul для композитного поля и для GF(2^8) по сути будут разные.

Матрицы изоморфного преобразования 11B => 13, N=1100
11B=>13 03 20 M =03A85C6870D2ACA0 M-1 =B1B042829AD45E54
11B=>13 03 22 M =737AF0C870D2ACA0 M-1 =A1B062A23A94BE14
11B=>13 03 35 M =91BC4AC4AE720CA0 M-1 =5FD0BAFAC21C2E9C
11B=>13 03 36 M =3FCE4664AE720CA0 M-1 =AFD01A5AE2DCCE5C
11B=>13 03 41 M =5DFC1C18DCAC72A0 M-1 =83700CACA4128692
11B=>13 03 45 M =81506EB8DCAC72A0 M-1 =B370CC6CE432E6B2
11B=>13 03 5B M =05E608CAA20CD2A0 M-1 =259024044C2A36AA
11B=>13 03 5E M =A7EADA6AA20CD2A0 M-1 =759064448C8A560A

Полученный результат позволяет выполнять операции композитного поля на инструкциях GFNI при изоморфном преобразовании из композитного поля в поле AES. И наоборот, позволяет перейти от операций в поле GF(2^8) к операциям меньшей разрядности GF(2^4). Такая техника может давать прирост производительности при реализации алгоритмов шифрования аппаратно и в микроконтроллерах.

В качестве тестового примера для проверки работоспособности методики выполняем генерацию S-Box AES в композитном поле с использованием изоморфного преобразования. При этом операция инверсии выполняется в арфиметике композитного поля.

S(x) = (A·M-1)·(M·x)-1 ⊕ C
A =F1E3C78F1F3E7CF816 - матрица аффинных преобразовний AES. 
C =6316
M =A7EADA6AA20CD2A016, -- матрица изоморфного преобразования
M-1 =759064448C8A560A16, -- обратное изоморфное преобразование
(A·M-1) =2F33DDCF49B6701E16

sbox[x] = affine(m_mul(A, M_), composite_inverse(affine(M, x)))^0x63;

Результат совпадает!

Расчет S-Box CLEFIA в композитных полях

Ранее мы рассматривали пример расчета нелинейного преобразования блочного шифра CLEFIA в изоморфных полях GF(2^8). Хотим закрепить навык выполнения расчетов в композитных полях и использовать технику преобразования в композитное поле, с расчетом инверсии на полиноме GF(2^4)/0x13.

Расчет нелинейного преобразования (биекции) дается следующим выражением:

S(x) = A2·(A1·x ⊕ C1)-1 ⊕ C2
A2,A1 - матрицы аффинных преобразовний.
инверсия считается в поле GF(2^8)/0x11D.
A1 =0x81605C6503015118 C1=0x1E
A2 =0x449002302058410A C2=0x69

Перейдем к расчетам в композитном поле:

S(x) = (A2·M-1)·((M·A1)·x ⊕ (M·C1))-1 ⊕ C2
инверсия считается в композитном поле GF((2^4)^2)/0x13.
M =B9C292428A441A2016
M-1 =9F547C4E5A805C0A16
(M·A1)  =FE297B311D0D060116,  (M·C1) =0x6A
(A2·M-1) =205054DA8048C31A16

// Расчет SBox - фрагмент кода
sbox[x] = affine(m_mul(A2,M_), composite_inverse(affine(m_mul(M,A1), x) 
                                               ^ affine(M,C1))) ^ C2;

Результат совпадает!

[1] D. Canright. A Very Compact S-box for AES
[2] A. Rudra et al. Efficient Rijndael Encryption Implementation with Composite Field Arithmetic

GF2p8 Аффинные преобразования (ч. 7)

Произведение матриц 8x8 бит

Инструкция аффинных преобразований gf2p8affine -- векторная. Аргументом является вектор элементов 64 бит, состоящий из матриц и вектор состоящий из октетов 8x8бит, над которыми выполняется преобразование. Мы хотим использовать эту инструкцию для перемножения матриц 8x8. Произведение матриц позволяет вычислять матрицы аффинного преобразования. В прошлом сообщении было показано, как выполнять транспонирование матриц. Операция аффинного преобразования обладает некоторой симметрией аргументов, которая позволяет менять аргументы местами при условии транспонирования.

gf2p8affine(A,B) = A·BT
gf2p8affine(A,BT) = A·B
gf2p8affine(AT,B) = (B·A)T= AT·BT

Чтобы инструкция gf2p8affine выдавала произведение матриц необходимо транспонировать второй аргумент.

// Произведение матриц A,B 8x8 бит
__m128i m_mul(__m128i A, __m128i B) {
	const __m128i E = _mm_set1_epi64x(0x0102040810204080);
	B = _mm_gf2p8affine_epi64_epi8(E, B, 0);// транспонирование
	A = _mm_gf2p8affine_epi64_epi8(A, B, 0);
	return A;
}

Результатом работы функции - произведение матриц A·B. Каждый аргумент состоит из двух матриц 8x8 по 64 бит, тип uint64x2_t. Параллельно можно выполнить два умножения матриц, независимо. Результат содержит две матрицы 8x8. Аналогично, при использовании инструкции 512 бит (_mm512_gf2p8affine_epi64_epi8), в результате получаем вектор uint64x8_t -- восемь матриц.

GF2p8 Аффинные преобразования (ч. 6)

Продолжаю тему, понял что некоторые разделы не раскрыты..

Умножение на константу в поле GF(2^8)

Нужен шаблон для матрицы умножения на константу в поле. Причем, умножение на константу реализовать надо в виде матрицы аффинного преобразования, чтобы иметь возможность последовательно накладывать и перемножать матрицы.

Матрица умножения состоит из столбцов, в которых располагаются степени числа 2. При умножении констант, должно быть выполнено редуцирование по образующему полиному.

M = [δ·27,δ·26,δ·25,δ·24,δ·23,δ·22,δ·21,δ·20]T

Матрицу нужно транспонировать, если заполнение производится по строкам.

// Алгоритм формирования матрицы умножения на константу в поле GF(2^8)
uint64_t m_mulC(uint8_t c, uint32_t Poly){
	uint64_t m = c;
	int i;
	for (i=1;i<8;i++){
		uint8_t v = gf2p8_mul(c, (1<<i), Poly);
		m = (m<<8) ^ v;
	}
	return m_transpose(m);
}

Редуцирование в поле GF(2^8)

Надо уметь выполнять редуцирование в форме аффинных преобразований - в одно действие. Дана старшая часть умножения без переноса, требуется выполнить редуцирование, т.е. перевести разряды в младщую часть. Редуцирование можно выполнить методом умножения на х8 mod Poly. В данном контексте матрица аффинного преобразования получается из выражения m_mulC(Poly & 0xFF, Poly), редуцирование - это матрица умножения в поле на специальную константу.

Таблица матриц редуцирования
0x11B: M=0xB1D3A6FD4B962C58
0x11D: M=0x71E2B51B478E1C38
0x165: M=0xDDBAA952A495F7EE
0x171: M=0x8D1A34685D37E3C6
0x177: M=0x4DD7E3C6C1CFD3A6
0x1C3: M=0x5BEDDAB468D0FBAD
0x1F5: M=0x2346AF5E9F1D1911

Транспонирование матрицы 8х8 бит в поле GF(2^8)

Трюк такой: можно поменять местами матрицу и вектора над которыми матрица работает, результат - можно получить транспонирование матрицы, если на входе диагональная.

// Транспонирование матриц 8x8 бит
__m128i transpose_8x8 (__m128i m){
    __m128i E = _mm_set1_epi64x(0x0102040810204080ULL);
    return _mm_gf2p8affine_epi64_epi8(E, m, 0);
}

Транспонирование матрицы M = 0x0123456789ABCDEF;
    M             MT
11101111 =EF   11110000 =F0
11001101 =CD   11001100 =CC
10101011 =AB   10101010 =AA
10001001 =89   00000000 =00
01100111 =67   11110000 =F0
01000101 =45   11001100 =CC
00100011 =23   10101010 =AA
00000001 =01   11111111 =FF

Чтобы убедиться, что это действительно транспонирование, применяем повторно. Повторное транспонирование возвращает матрицу к исходному виду. Транспонирование в таком виде нужно, чтобы динамически составлять матрицы умножения на константу. Ниже привожу функцию динамической генератции матрицы аффинного преобразования, которая выполняет умножение на константу в поле GF(2^8)/0x11B

// Составление Матрицы умножения на константу в поле 0x11B
__m128i m_mulC_11B(uint8_t c) {
    __m128i E = _mm_set1_epi64x(0x0102040810204080ULL);
    __m128i m = _mm_gf2p8mul_epi8(E, _mm_set1_epi8(c));
    return _mm_gf2p8affine_epi64_epi8(E, m, 0);// транспонирование 8x8
}

Матрица редуцирования [1B, 0x11B] A= 0xB1D3A6FD4B962C58
Матрица умножения на  [A5, 0x11B] A= 0xABFCF9581A356AD5

GF2p8 Аффинные преобразования (ч. 4)

Вычисление CRC на операции умножения в изоморфном представлении поля GF(2^8)

Мы приблизились к тому, ради чего все затевали.. Нет, вычисление CRC-8 не самоцель. Это самый простой пример применения арифметики Галуа, на котором можно проверить методику. Наша цель - изоморфные преобразования поля, чтобы свести все варианты расчета к одному полиному 0x11B.

Хотелось бы донести хотя бы до одного человека, как выводить алгоритмы. На первом шаге берем операцию, на которой работает CRC:
a(x)×x⁸ mod P(x) -- операция '×'-умножение без переноса или сдвиг.
Это можно представить иначе:
a(x)⊗(x⁸ mod P(x)) -- операция '⊗'-умножение в поле P(x).
Для полинома 0x11D выражение принимает вид: a(x)⊗0x1D -- операция '⊗'-умножение в поле P(x)=0x11D.

// CRC-8 в поле 0x11D
uint8_t crc = 0xFD;// начальное значение
uint8_t x8  = 0x1D;
for (i=0; i<length; i++) {
    crc = gf2p8_mul(crc ^ data[i], x8, 0x11D);
}
return crc;

M^-1·M(a(x)×0x1D mod 0x11D) = M^-1((M·a(x)) × (M·0x1D) mod 0x11B) -- операция '×'-сначала выполняется в поле P(x)=0x11D, а после преобразования в поле P(x)=0x11B.

// CRC-8 в изоморфном представлении 0x11D=>0x11B
// матрицы изоморфного преобразования {11D, 02} =>{11B, FB}
uint64_t M =0x57EED8E22A228202;
uint64_t Mr=0x4F803A301CA0E8C0;

uint8_t crc = affine_byte(M, 0xFD);// начальное значение
uint8_t x8  = affine_byte(M, 0x1D);// преобразование 0x11D=>0x11B
for (i=0; i<length; i++) {
    crc^= affine_byte(M, data[i]);
    crc = gf2p8_mul(crc, x8, 0x11B);
}
return affine_byte(Mr, crc);// обратное преобразование 0x11B=>0x11D

FD => 6A        57EED8E22A228202
1D => 11        57EED8E22A228202
DE => 7E        4F803A301CA0E8C0
CRC8I=7E ..ok

Операция вычисляется на двух инструкциях в цикле: gf2p8_mul и affine. Работает!!