CODEWITHSUNDEEP
algorithm
Zagorax
Zagorax

Reputation: 11900

Find the inverse of a number modulo a prime

I was thinking about an algorithm to solve the congruence ax = 1 mod p with p prime. I was thinking about using Fermat theorem. Since I know that

a ^ (p-1) = 1 mod p

and that

a ^ (p-1) = a * (a ^ (p-2))

It means that a ^ (p-2) mod p is the solution. Unfortunately this solution, although mathematically correct, isn't good for computer since for big primes I have to do a ^ (p-2) which is usually not calculable.

Which algorithm is good for computer science?

Upvotes: 3

Views: 11588

Answers (3)

IVlad
IVlad

Reputation: 43517

There is no reason this isn't a good algorithm for computers, you just need to be careful about the implementation, which isn't exactly trivial I guess, but it's not difficult either.

Just use exponentiation by squaring, then it most likely won't matter how big p is.

a^n = a^(n / 2) * a^(n / 2) for n even
    = a*a^(n - 1) for n odd

Upvotes: 1

Daniel Fischer
Daniel Fischer

Reputation: 183978

since for big primes I have to do a ^ (p-2) which is usually not calculable.

You need modular exponentiation, so with the exponentiation by squaring mentioned by IVlad you only need Θ(log p) modular multiplications of numbers of size at most p-1. The intermediate results are bounded by p^2, so despite a^(p-2) not being calculable for large primes, (a ^ (p-2)) % p usually is. That method is simple to implement:

unsigned long long invert_mod(unsigned long long a, unsigned long long p) {
    unsigned long long ex = p-2, result = 1;
    while (ex > 0) {
        if (ex % 2 == 1) {
            result = (result*a) % p;
        }
        a = (a*a) % p;
        ex /= 2;
    }
    return result;
}

but has a few drawbacks. (p-1)^2 must be representable in the used type (not a problem [except for huge p] if arbitrary precision integers are used), or you get invalid results due to overflow, and it always uses at least log (p-2)/log 2 modular multiplications.

The extended Euclidean algorithm, as suggested by user448810, or equivalently the continued fraction method, never produces intermediate values larger than p, thus avoiding all overflow problems if p is representable, and usually needs fewer divisions. Additionally, it computes the modular inverse not only for primes, but for any two coprime numbers.

unsigned long long invert_mod(unsigned long long a, unsigned long long p) {
    unsigned long long new = 1, old = 0, q = p, r, h;
    int pos = 0;
    while (a > 0) {
        r = q%a;
        q = q/a;
        h = q*new + old;
        old = new;
        new = h;
        q = a;
        a = r;
        pos = !pos;
    }
    return pos ? old : (p - old);
}

The code is slightly longer, but an optimising compiler ought to compile it to a short loop using only one division per iteration.

Upvotes: 14

user448810
user448810

Reputation: 17866

The normal way to compute the modular inverse is by the extended Euclidean algorithm:

function inverse(x, m)
    a, b, u := 0, m, 1
    while x > 0
        q, r := divide(b, x)
        x, a, b, u := b % x, u, x, a - q * u
    if b == 1 return a % m
    error "must be coprime"

Upvotes: 8

Related Questions