Calculate * mod without overflow

Question

I need to compute a*a mod n but a quite large, which leads to overflow when I square it. Executing ((a % n)*(a % n)) % n does not work because (n-1) ² may overflow. This is in C ++ and I use int64_t

Edit:

Example values: a = 821037907258 and n = 800000000000, which overflows if you square it.

I am using DevCPP and I have already tried to get integer libraries to work to no avail.

Edit 2:

No, there are no patterns for these numbers.

Answer 1

If you cannot use a library with a large integer, and you do not have a native uint128_t (or similar), you will need to do this manually.

One option is to express a as the sum of two 32-bit values, i.e. a = 2 ³² b + c, where b contains 32 msbs and c contains 32 lbs. Quadraging is a set of four cross-multiplications; each result is guaranteed to fit into the 64-bit type. Then you do the modulo operation when you recombine individual terms (carefully taking into account the shifts needed to rebuild everything).

Answer 2

I know that you no longer need this, and there is an alternative solution, but I want to add an alternative method for its implementation. It offers two different methods: double and add an algorithm and a mod(a + b, n) processing method with overflow detection.

The double addition algorithm is usually used in fields where multiplication is impossible or too expensive to calculate directly (e.g. elliptic curves), but we could take it to handle it in our situation, to handle overflows instead.

The following code is probably slower than the decision made (even when optimizing it), but if the speed is not critical, you may prefer it for clarity.

 unsigned addmod(unsigned x, unsigned y, unsigned n) { // Precondition: x<n, y<n // If it will overflow, use alternative calculation if (x + y <= x) x = x - (n - y); else x = (x + y) % n; return x; } unsigned sqrmod(unsigned a, unsigned n) { unsigned b; unsigned sum = 0; // Make sure original number is less than n a = a % n; // Use double and add algorithm to calculate a*a mod n for (b = a; b != 0; b >>= 1) { if (b & 1) { sum = addmod(sum, a, n); } a = addmod(a, a, n); } return sum; } 

Answer 3

Here is an implementation with the double addition of a * b % m multiplication, without overflow, regardless of the size of a, b and m.

(Note that the lines res -= m and temp_b -= m rely on a 64-bit unsigned integer overflow to produce the expected results. This should be good, because an unrecognized integer overflow is correctly defined in C and C ++. The reason is important for use unsigned integers .)

 uint64_t mulmod(uint64_t a, uint64_t b, uint64_t m) { uint64_t res = 0; uint64_t temp_b; /* Only needed if b may be >= m */ if (b >= m) { if (m > UINT64_MAX / 2u) b -= m; else b %= m; } while (a != 0) { if (a & 1) { /* Add b to res, modulo m, without overflow */ if (b >= m - res) /* Equiv to if (res + b >= m), without overflow */ res -= m; res += b; } a >>= 1; /* Double b, modulo m */ temp_b = b; if (b >= m - b) /* Equiv to if (2 * b >= m), without overflow */ temp_b -= m; b += temp_b; } return res; } 

This is my modification of another answer to another similar question .

Answer 4

No need for a math library of arbitrary precision. This is easy to do without double-width arithmetic and without loops

• If a <2 ^N-1, we can just multiply it directly without worrying about overflow
• If a = 2 ^N-1, then a ² 2 ^N 2 ^N - n + n 2 ^N - n (mod n)

• If a> 2 ^N-1, then let A̅ be two additions to a, that is, A̅ = -a = 2 ^N - a <2 ^N-1 . This means that A̅ ² will not overflow, and we can calculate ² % n using the expression below

  a ² ≡ [2 ^N - (2 ^N - a)] ² ≡ (2 ^N - a) ² + 2 ^2N - 2.2 ^N (2 ^N - a) ≡ A̅ ² + 2 ^2N - 2.2 ^N A̅ ( mod n)

  • To get 2 ^2N % n, we need to divide it into different cases:
    • If n <2 ^N-1, then 2 ^N × 2 ^N % n = (2 ^N % n) × (2 ^N % n)% n does not overflow
    • If n> 2 ^N-1, then 2 ^N - n <2 ^N-1 and 2 ^N % n = (2 ^N - n + n)% n = (2 ^N - n)% n does not overflow
    • If n = 2 ^N-1, then this is trivial
  • To get 2.2 ^N A̅% n, we will do the same. Please note that 2A̅ does not overflow
    • If n <2 ^N-1, then (2A̅% n) (2 ^N % n)% n does not overflow
    • If n> 2 ^N-1, then use (X - n + n)% n = (X - n)% n for both terms as above to avoid overflow
    • If n = 2, ^N-1 : trivial

The result is

 template<typename T> T squareMod(T a, T n) { constexpr auto N = sizeof(T)*CHAR_BIT; constexpr T midpoint = T(1) << N/2; if (a < midpoint) return a*a % n; // use the normal way, because a*a can't overflow else if (a == midpoint) return (std::numeric_limits<T>::max() - n + 1) % n; else { TA = -a; T x = A*A % n; // a is large, so -a is small and A*A won't overflow T y, z; if (n == midpoint) return x; // 2ᴺ % 2ᴺ⁻¹ = 0 else if (n < midpoint) { T twoPowN_modN = std::numeric_limits<T>::max() % n; if (twoPowN_modN == n) twoPowN_modN = 0; y = twoPowN_modN*twoPowN_modN % n; z = (2*A % n)*twoPowN_modN % n; } else { T twoPowN_modN = -n; y = twoPowN_modN*twoPowN_modN % n; z = ((2*A - n) % n)*twoPowN_modN % n; } return (x - z + y) % n; } } 

Answer 5

You can implement the multiplication yourself by adding n of each run and immediately changing the result.

Answer 6

I really believe that ((a% n) * (a% n))% n should work for positive integers. Why do you think this is not working? Do you have a counter example? If n can be negative, then the% undefined operator.