Naglasimo za početak da je u programu potrebno množiti dva broja koji mogu biti do \(10^9\) i izračunati proizvod koji može biti do \(10^{18}\). Zato je u programu potrebno koristiti 64-bitno zapisane brojeve. U jeziku C++ možemo koristiti tip long long. U jeziku C# možemo koristiti tip long.

Pretraga grubom silom

Za svaki učitani broj \(x\) treba odrediti \(y\) takav da je \(a\cdot_n y = x\). Broj \(y\) se, dakle, može odrediti tako što se izvrši deljenje broja \(x\) brojem \(a\) po modulu \(n\). Deljenje po modulu je veoma specifična operacija, drugačija od svih operacija po modulu koje smo do sada sreli. Čak ni postojanje količnika po modulu \(n\) nije uvek garantovano (ako su brojevi \(a\) i \(n\) uzajamno prosti, kao što jesu u ovom zadatku tada modularni količnik postoji za svako \(x\) i jedinstven među brojevima od \(0\) do \(n-1\)).

Modularno deljenje grubom silom podrazumevalo bi primenu algoritma pretrage za svako \(x\) tj. isprobavanje svih mogućih brojeva \(y\) od \(0\) do \(n-1\) i proveru da li je \(a\cdot_n y = x\). Pošto \(n\) može biti prilično veliki broj, pronalaženje inverza i za jedno \(x\) može biti veoma neefikasno.

#include <iostream>

using namespace std;

// pokusavamo da podelimo x sa a po modulu n
// ako postoji rezultat se smesta u y, a funkcija vraca podatak o tome
// da li kolicnik postoji
bool deli(long long x, long long a, long long n, long long& y) {
  // resenje grubom silom
  // posto je (a * y) mod n = (a mod n * y mod n) mod n
  // mozemo a zameniti sa a % n
  a %= n;
  // proveravamo sve moguce kandidate za y
  for (y = 0; y < n; y++)
    // ako je 
    if ((a * y) % n == x)
      return true;
  // kolicnik ne postoji
  return false;
}

int main() {
  // ucitavamo brojeve
  long long a, n;
  cin >> a >> n;
  // ucitavamo sve brojeve do kraja ulaza
  long long x;
  while (cin >> x) {
    // delimo ih sa a po modulu n (po uslovima zadatka kolicnik
    // postoji)
    long long y;
    deli(x, a, n, y);
    cout << y << endl;
  }
  return 0;
}

Umesto nezavisnog deljenja svakog broja \(x\) brojem \(a\) (po modulu \(n\)), bolji pristup je izračunavanje modularnog inverza broja \(a\) i kasnije množenje (po modulu \(n\)) svih brojeva \(x\) tim inverzom

Za dati broj \(a\) potrebno je pronaći njegov inverz po modulu \(n\) tj. broj \(a^{-1}\) takav da je \(a^{-1} \cdot_n a = 1\), tj. da je \((a^{-1} \cdot a)\,\mathrm{mod}\,n = 1\). Ako se takav broj nađe, onda važi da je \(a^{-1} \cdot_n (a \cdot_n y) = a^{-1}\cdot_n x\). Pošto je \(a^{-1} \cdot_n (a \cdot_n y) = (a^{-1}\cdot a) \cdot_n y = 1 \cdot_n y = y\), važi da je \(y = a^{-1}\cdot_n x\), što nam omogućava da dešifrovanje tj.  deljenje po modulu izvršimo operacijom množenja po modulu.

Određivanje modularnog inverza se može uraditi grubom silom tako što se redom ispitaju svi brojevi \(a^{-1}\) od \(0\) do \(n-1\) i proveri se da li je \(a^{-1}\cdot_n a = 1\). Iako se pretraga ovaj put sprovodi samo jednom, a ne za svako \(x\) iznova, ovo rešenje je opet veoma neefikasno u slučaju kada je \(n\) veliki broj.

#include <iostream>

using namespace std;

// odredjujemo inverz broja a po modulu n i smestamo ga u promenljivu t
// funkcija vraca podatak o tome da li modularni inverz postoji
bool inverz(long long a, long long n, long long& t) {
  // resenje grubom silom

  // posto je (a * t) mod n = (a mod n * t mod n) mod n
  // mozemo a zameniti sa a % n
  a %= n;
  // proveravamo sve moguce kandidate za inverz
  for (t = 0; t < n; t++)
    // ako je ostatak 1, nasli smo inverz
    if ((t * a) % n == 1)
      return true;
  // inverz ne postoji
  return false;
}

int main() {
  // ucitavamo brojeve
  long long a, n;
  cin >> a >> n;
  // trazimo inverz elementa a po modulu n (po uslovima zadatka on
  // mora da postoji)
  long long inv_a;
  inverz(a, n, inv_a);
  // ucitavamo sve brojeve do kraja ulaza
  long long x;
  while (cin >> x)
    // delimo ih sa a po modulu tako sto ih mnozimo sa inverzom od a
    cout <<  (inv_a * x) % n << endl;
  return 0;
}

Prošireni Euklidov algoritam i Bezuova teorema

Jedna izrazito važna posledica Euklidovog algoritma opisanog u zadatku [euklid] je Bezuova teorema koja kaže da za svaka dva prirodna broja \(a\) i \(b\) postoje celi brojevi \(s\) i \(t\) tako da važi da je \(s\cdot a + t \cdot b = nzd(a, b)\). Brojevi \(s\) i \(t\) se mogu izračunati tzv. proširenim Euklidovim algoritmom. Posmatrajmo narednu tabelu.

\[ \begin{array}{rclrclrcl} r_0 &=& a & s_0 &=& 1& t_0 &=& 0\\ r_1 &=& b & s_1 &=& 0& t_1 &=& 1\\ r_2 &=& r_0 - q_1r_1 &s_2 &=& s_0 - q_1s_1 &t_2 &=& t_0 - q_1t_1\\ \ldots\\ r_{i+1} &=& r_{i-1} - q_ir_i & s_{i+1} &=& s_{i-1} - q_is_i & t_{i+1} &=& t_{i-1} - q_it_i\\ \ldots \\ r_{k} &=& r_{k-2} - q_{k-1}r_{k-1} & s_{k} &=& s_{k-2} - q_{k-1}s_{k-1} & t_{k} &=& t_{k-2} - q_{k-1}t_{k-1}\\ r_{k+1} &=& 0 & s_{k+1} &=& s_{k-1} - q_{k}s_{k} & t_{k+1} &=& t_{k-1} - q_{k}t_{k} \end{array} \]

U svakom njenom koraku važi da je \(s_i\cdot a + t_i\cdot b = r_i\). Zaista, u prva dva koraka jednakost je trivijalno ispunjena. Ako pretpostavimo da važi \(s_{i-1}a + t_{i-1}b = r_{i-1}\) i \(s_ia + t_ib = r_i\), tada je na osnovu tabele \(r_{i+1} = (s_{i-1}a + t_{i-1}b) - q_i(s_ia + t_ib) = (s_{i-1} - q_is_i)\cdot a + (t_{i-1}-q_it_i)\cdot b\). Pošto je u tabeli stavljeno da je \(s_{i+1} = s_{i-1} - q_is_i\) i \(t_{i+1} = t_{i-1} - q_it_i\), i u koraku \(i+1\) važi da je \(s_{i+1}a + t_{i+1}b = r_{i+1}\). Dakle, važi da je \(s_{k}a + t_{k}b = r_{k}\), a pošto je \(r_k\) NZD brojeva \(a\) i \(b\), brojevi \(s_k\) i \(t_k\) su traženi Bezuovi koeficijenti. Činjenica da je \(r_k\) NZD brojeva \(a\) i \(b\) dokazana je u zadatku [euklid].

Dodatno, važi da je \(s_{k+1} = \pm \frac{b}{nzd(a, b)}\) i \(t_{k+1} = \pm \frac{a}{nzd(a, b)}\). Zaista, važi da su brojevi \(s_{k+1}\) i \(t_{k+1}\) uzajamno prosti (to sledi iz identiteta \(s_kt_{k+1} - t_{k}s_{k+1} = (-1)^k\)). Pošto uz to važi i da je \(as_{k+1} + bt_{k+1} = 0\), važi da broj \(s_{k+1}\) deli broj \(b\) dok \(t_{k+1}\) deli broj \(a\), pa tvrđenje sledi. Može se pokazati i da brojevi \(s_i\) i \(t_i\) obrazuju serije alternirajućeg znaka koje su (posle nekih prvih članova) strogo rastuće po apsolutnoj vrednosti. Zato važi da je \(|s_i| < \frac{b}{nzd(a, b)}\) i \(|t_i| < \frac{a}{nzd(a, b)}\), što je važno svojstvo koje nas osigurava od nastanka prekoračenja tokom izračunavanja.

Implementacija je jednostavno uopštenje osnovne implementacije Euklidovog algoritma. U po dve promenljive održavamo dve uzastopne vrednosti niza \(r_i\), niza \(s_i\) i niza \(t_i\), u petlji ih ažuriramo na osnovu prikazanih formula sve dok član niza \(r_i\) ne postane jednak nuli.

Videćemo da se prošireni Euklidov algoritam može primeniti u mnogo različitih situacija.

Vratimo se našem originalnom problemu. Zadatak smo sveli na određivanje modularnog inverza broja \(a\) po modulu \(n\) tj. na to da treba pronaći broj \(q\) tako da je \(a^{-1} \cdot a = q\cdot n + 1\), tj. \(-q\cdot n + a^{-1} \cdot a = 1\). Pošto su brojevi \(a\) i \(n\) uzajamno prosti po pretpostavci zadatka (u suprotnom ne bismo mogli da garantujemo postojanje rešenja), važi da je NZD brojeva \(a\) i \(n\) jednak 1. Na osnovu Bezuove teoreme moguće je pronaći brojeve \(s\) i \(t\) takve da je \(s\cdot n + t \cdot a = 1\). Kada proširenim Euklidovim algoritmom pronađemo te keoficijente, koeficijent \(t\) je traženi modularni inverz. Ako je taj broj negativan, onda se za inverz može uzeti \(t + n\) (jer tada važi \((s-a)\cdot n + (t+n)\cdot a = 1\), a pošto je \(|t| < \frac{n}{nzd(a, n)} = n\), broj \(t+n\) je sigurno pozitivan).

#include <iostream>

using namespace std;

long long bezu(long long a, long long b, long long& s, long long& t) {
  long long r0 = a, r1 = b;
  long long s0 = 1, s1 = 0;
  long long t0 = 0, t1 = 1;
  while (r1 > 0) {
    long long q = r0 / r1;
    long long tmp;
    
    tmp = r0;
    r0 = r1;
    r1 = tmp - q*r1;
    
    tmp = s0;
    s0 = s1;
    s1 = tmp - q*s1;
    
    tmp = t0;
    t0 = t1;
    t1 = tmp - q*t1;
  }
  s = s0; t = t0;
  return r0;
}

bool inverz(long long a, long long n, long long& t) {
  long long s;
  long long r = bezu(n, a, s, t);
  if (t < 0) t += n;
  return r == 1;
}

int main() {
  long long a, n;
  cin >> a >> n;
  long long inv_a;
  inverz(a, n, inv_a);
  long long x;
  while (cin >> x)
    cout <<  (inv_a * x) % n << endl;
  return 0;
}

Primetimo da nam vrednost broja \(s\) zapravo nije potrebna, pa se njeno izračunavanje u ovoj specijalnoj primeni proširenog Euklidovog algoritma zapravo može izostaviti.

#include <iostream>

using namespace std;

// izracunavamo broj t takav da je a * t mod n = 1
// vracamo podatak o tome da li takav broj t postoji
bool inverz(long long a, long long n, long long& t) {
  long long r = n, r1 = a;
  t = 0; long long t1 = 1;
  while (r1 > 0) {
    long long q = r / r1;
    long long tmp;
    
    tmp = r;
    r = r1;
    r1 = tmp - q*r1;
    
    tmp = t;
    t = t1;
    t1 = tmp - q*t1;
  }
  if (t < 0) t += n;
  return r == 1;
}


int main() {
  long long a, n;
  cin >> a >> n;
  long long inv_a;
  inverz(a, n, inv_a);
  long long x;
  while (cin >> x)
    cout <<  (inv_a * x) % n << endl;
  return 0;
}

Ojlerova teorema

Još jedan način da efikasno izračunamo modularni inverz je korišćenje Ojlerove teoreme. Neka \(\varphi(n)\) označava broj brojeva koji su uzajamno prosti sa brojem \(n\). Ojlerova funkcija način njenog efikasnog izračunavanja opisan je u zadatku [ojlerova_funkcija]. Ojlerova teorema tvrdi da ako su brojevi \(a\) i \(n\) uzajmno prosti važi da je \(a^{\varphi(n)} \,\mathrm{mod}\,n = 1\) tj. da je \(\left(a^{\varphi(n)-1}\right)\cdot_n a \,\mathrm{mod}\,n = 1\). Zato je \(a^{-1}\) po modulu \(n\) jednak broju \(a^{\varphi(n)-1}\). Pošto eksponent može biti veliki broj stepenovanje treba izvršiti nekim efikasnim algoritmom stepenovanja (stepenovanje se, naravno, vrši po modulu \(n\)).

#include <iostream>

using namespace std;

// broj brojeva manjih od n koji su uzajamno prosti sa n
long long phi(long long n) {
  // Ojlerova formula kaze da je taj broj jednak
  // n * proizvod(1 - 1/p) za sve proste brojeve p koji dele broj n
  long long d = 2;
  long long proizvod = n;
  // dok broj n ne postane prost
  while (d*d <= n) {
    if (n % d == 0) {
      // nasli smo nov prost faktor pa azuriramo proizvod
      proizvod = (proizvod / d) * (d - 1);
      // uklanjamo ostala pojavljivanja tog faktora
      while (n % d == 0)
        n /= d;
    }
    // prelazimo na sledeceg kandidata
    d++;
  }
  // n je prost broj, pa treba azurirati proizvod i za njega
  if (n > 1)
    proizvod = (proizvod / n) * (n - 1);
  // vracamo rezultat
  return proizvod;
}

// izracunavamo x^k mod n koristeci algoritam efikasnog stepenovanja
long long stepen_po_modulu(long long x, long long k, long long n) {
  k %= n;
  long long stepen = 1;
  while (k > 0) {
    if (k % 2 == 1)
      stepen = (stepen * x) % n;
    x = (x * x) % n;
    k /= 2;
  }
  return stepen;
}

// izracunavamo broj t takav da je a * t mod n = 1
// funkcija pretpostavlja da su a i n uzajamno prosti i tada t uvek postoji
long long inverz(long long a, long long n) {
  return stepen_po_modulu(a, phi(n)-1, n);
}


int main() {
  // ucitavamo brojeve
  long long a, n;
  cin >> a >> n;
  // trazimo inverz elementa a po modulu n (po uslovima zadatka on
  // mora da postoji)
  long long inv_a = inverz(a, n);
  // ucitavamo sve brojeve do kraja ulaza
  long long x;
  while (cin >> x)
    // delimo ih sa a po modulu tako sto ih mnozimo sa inverzom od a
    cout <<  (inv_a * x) % n << endl;
  return 0;
}