Линейное преобразование в ГОСТ Р 34.12-2015

Здравствуйте!
Мне хочется написать реализацию недавно принятого на вооружение симметричного блочного алгоритма шифрования "Кузнечик". Сейчас изучаю стандарт, пытаюсь посчитать каждый шаг алгоритма "на пальцах" чтобы подтвердить своё понимание, и один момент никак не сходится с тестовым примером.
Текст ГОСТа представлен на сайте ТК 26: http://tc26.ru/standard/gost/GOST_R_3412-2015.pdf
О чём речь?
Проблема в трансформации R, которая преобразует 16-байтный массив, а точнее — в l, которая из этого массива получает всего один байт данных.
R(a[]) = l(a[]) & a[15] & a[14] & .. & a[1] l(a[]) = 148 * a[15] + 32 * a[14] + ... + 1 * a[0]
Операции сложения и умножения осуществляются в поле F, а константы являются элементами поля в указанном ранее смысле.
Как я понял, речь идёт о сложении и умножении по модулю 255
F конечное поле GF(2)[x]∕p(x), где p(x) = x8 + x7 + x6 + x + 1 ∈ GF(2)[x]
Что случилось?
В приложении к стандарту приводятся контрольные примеры для преобразований, используемых в алгоритме шифрования, и для R в том числе.
R(00000000000000000000000000000100) = 94000000000000000000000000000001 R(94000000000000000000000000000001) = a5940000000000000000000000000000
Шестнадцатеричные числа — представление массива из 16 байт, причём первый (а точнее, нулевой) элемент находится в конце строки.
Итак, R(a[]) = l(a[]) & a[15] & ... & a[1]
В первом примере всё хорошо:
l(a[]) = 148 * a[1] = 148 = 0x94 R(a[]) = 0x94 & a[15..1] = 0x94, 0x00, ... 0x01
Во втором случае всё сложнее. Если я правильно понял про сложение и умножение, озвученное выше, то:
l(a[]) = 148 * a[15] + 1 * a[0] = (148 * 148) mod 255 + 1 * 1 = 229 + 1 = 230 = E6
R(a[]) = 0xE6, 0x94, 0x00, ...
Стоп! В контрольном примере старшие байты получились 0xA5, 0x94!
Что делать?
Собственно, вот и сам вопрос :)
В какое место в моём потоке сознания прокралась ошибка и результат функции во l() втором примере оказался не тем, каким надо, и что сделать, чтобы всё закончилось хорошо?
Обновление №1
Как оказалось, умножение в поле — совсем не "просто умножить и взять остаток от деления", а целая процедура, подробно описанная в GF(256) finite field multiplication function in C#
Более того, сложение в поле GF(2^n) оказалось простым XOR'ом.
Итого для второго примера получилось:
l(a[]) = (148 * 0x94) + (1 * 0x01) = 0x91 + 0x01 = 0x90
Что в общем-то всё ещё далеко от A5
Привожу код расчёта l()
byte[] l_mul = { 1, 148, 32, 133, 16, 194, 192, 1, 251, 1, 192, 194, 16, 133, 32, 148 };
public byte l (byte[] a) { byte p = 0x00; for (int i = 0; i < 16; i++) { p ^= gf.Mul(l_mul[i], a[i]); } return p; }

Ответ

В общем, проблема была в процедуре умножения в поле. Если результат умножения полиномов получается слишком большой и не входит в GF(256), то его надо поделить на некоторый неприводимый полином.
В примере из Википедии и упомянутом ранее вопросе на Stack Overflow используется полином x^8 + x^4 + x^3 + x + 1 (100011011b = 0x11b или 0x1b), который применяется в AES.
В случае же ГОСТа из определения поля F необходимо использовать другой полином: x^8 + x^7 + x^6 + x + 1 = 0x1c3 или 0xc3. Теперь с контрольными примерами всё сходится :)
public byte MulSlow (byte a, byte b) { byte p = 0; byte counter; byte hi_bit_set; for (counter = 0; counter < 8 && a != 0 && b != 0; counter++) { if ((b & 1) != 0) p ^= a; hi_bit_set = (byte)(a & 0x80); a <<= 1; if (hi_bit_set != 0) a ^= 0xc3; /* x^8 + x^7 + x^6 + x + 1 */ b >>= 1; } return p; }