Быстрое умножение по модулю 2^16 + 1

Question

Шифр ​​IDEA использует умножение по модулю 2^16 + 1. Есть ли алгоритм для выполнения этой операции без общего оператора modulo (только по модулю 2^16 (усечение))? В контексте IDEA нуль интерпретируется как 2^16 (это означает, что ноль не является аргументом нашего умножения, и это не может быть результатом, поэтому мы можем сохранить один бит и сохранить значение 2^16 как битовый шаблон 0000000000000000). Мне интересно, как эффективно его реализовать (или вообще возможно) без использования стандартного оператора modulo.

Answer 1

Вы можете использовать тот факт, что (N-1)% N == -1.

Таким образом, (65536 * а)% 65537 == -a% 65537.
Кроме того, -a% 65537 == -a + 1 (по модулю 65536), при 0 интерпретируется как 65536

uint16_t fastmod65537(uint16_t a, uint16_t b) 
{ 
    uint32_t c; 
    uint16_t hi, lo; 
    if (a == 0) 
     return -b + 1; 
    if (b == 0) 
     return -a + 1; 

    c = (uint32_t)a * (uint32_t)b; 
    hi = c >> 16; 
    lo = c; 
    if (lo > hi) 
     return lo-hi; 
    return lo-hi+1; 
} 

Единственная проблема здесь, если hi == lo, то результат будет 0. К счастью набор тестов подтверждает, что он на самом деле не может быть ...

int main() 
{ 
    uint64_t a, b; 
    for (a = 1; a <= 65536; a++) 
     for (b = 1; b <= 65536; b++) 
     { 
      uint64_t c = a*b; 
      uint32_t d = (c % 65537) & 65535; 
      uint32_t e = m(a & 65535, b & 65535); 
      if (d != e) 
       printf("a * b % 65537 != m(%d, %d) real=%d m()=%d\n", 
         (uint32_t)a, (uint32_t)b, d, e); 
     } 
    } 

Выход: нет

Answer 2

Во-первых, случай, когда a или b равен нулю. В этом случае, он интерпретируется как имеющий значение 2^16, поэтому элементарно по модулю арифметика говорит нам, что:

result = -a - b + 1; 

, потому что (в контексте IDEA) мультипликативной инверсией 2^16 еще 2^16, а его младшие 16 бит - все нули.

Общий случай гораздо проще, чем кажется, теперь, когда мы взяли на себя "0" частный случай (2^16 + 1 является 0x10001):

/* This operation can overflow: */ 
unsigned result = (product & 0xFFFF) - (product >> 16); 
/* ..so account for cases of overflow: */ 
result -= result >> 16; 

Собираем вместе:

/* All types must be sufficiently wide unsigned, e.g. uint32_t: */ 
unsigned long long product = a * b; 
if (product == 0) { 
    return -a - b + 1; 
} else { 
    result = (product & 0xFFFF) - (product >> 16); 
    result -= result >> 16; 
    return result & 0xFFFF; 
}