What algorithm is prefered to do a x b mod P with big numbers (256 bits)

Question

I'm trying to implement multiple precision arithmetic operations modulo P, with P < 2^256. More specifically, P = 2^256 - 2^32 - 977.

I want to support the following operations: +, -, *, /, pow (each mod P)

As P is close to 2^256, numbers are represented with 8 u32 or 4 u64.

a + b mod P can be done like this (in pseudo code):

n = a + b
if overflow: # i.e. over 2^256
    # add 2^256 - P to come back modulo P
    n += 2**32 + 977
else:
    if n >= P:
        # P <= n <= 2^256
        n -= P

--

For a * b mod P, my first intention was to simply do a long multiplication but that seems slow as I would need the carry to be 256 bits as well.

Are there any recommended algorithms to calculate a * b modulo P efficiently (using arrays of u32 / u64)?

I'm mostly interested in the multiplication because:

• a^x mod P can be an optimized version of a * a * ... * a mod P
• a / b mod P can be calculated as a * b^{P-2} using fermats little theorem

Note: Bitcoin implements these operations with numbers represented with 10 x uint26 instead of 8 uint32 so each "digit" keeps 6 bits but I'm not familiar with their methods.

Answer 1

Here is one reasonable method:

To multiply a 32-bit integer by a 256-bit integer modulo $P$, multiply the integers using arbitrary-precision arithmetic (see How do computers perform operations on numbers that are larger than 64 bits? this can be done with 8 32x32 -> 64 multiplications, and then some 32-bit additions) to get a 288-bit product, then reduce the product modulo $P$ (divide the product by $P$, and keep the remainder).

To multiply a 256-bit integer $X$ by a 256-bit integer $Y$ modulo $P$, write

$$X = 2^{224} X_7 + \dots + 2^{32} X_1 + X_0,$$

then do the following:

• set $B := 0$
• for $i := 7,6,\dots,0$:
  • set $A := X_i \times Y \bmod P$
  • set $B := 2^{32} \times B \bmod P$
  • set $B := A + B \bmod P$

At the end, $B$ will hold the product $A \times B \bmod P$. Each step can be computed using the method of the first paragraph of this answer as it only involves 32x256 -> 288 modular multiplications.

Answer 2

Note that with your particular P, $a \cdot 2^{256} + b \mod P = $ $ b + a \cdot (2^{32}+977) \mod P$. The result will only rarely be P or slightly larger, in which case you subtract P once more.