Atenţie! Aceasta este o versiune veche a paginii, scrisă la 2008-01-31 15:24:07.
Revizia anterioară Revizia următoare

Lucrul cu numere mari

De multe ori, in probleme, apar situatii cand este nevoie sa memoram numere intregi foarte mari (de ordinul sutelor de cifre), iar uneori trebuie sa efectuam si operatii aritmetice cu aceste numere.

Numerele reprezentate pe vectori le vom numi pur si simplu "vectori", iar numerele reprezentate printr-un tip ordinal de date le vom numi, printr-o analogie usor fortata cu matematica, "scalari". Sa vedem acum cum se efectueaza operatiile elementare pe aceste numere.

Initializarea

Un vector poate fi initializat in trei feluri: cu 0, cu un scalar sau cu un alt vector.

La initializarea cu 0, singurul lucru pe care il avem de facut este sa setam numarul de cifre pe 0. De aceea, este practic inutil sa implementam aceasta functie ca atare; putem folosi in loc singura instructiune pe care ea o contine.

void Atrib0(Huge H) {
  H[0] = 0;
}

La initializarea cu un scalar nenul, trebuie sa asezam fiecare cifra pe pozitia corespunzatoare, afland in paralel si numarul de cifre. Se incepe cu cifra unitatilor, si la fiecare pas se pune in vector cifra cea mai putin semnificativa, dupa care numarul de reprezentat se imparte la 10 (neglijandu-se restul), iar numarul de cifre se incrementeaza.

void AtribValue(Huge H, unsigned long X) {
  H[0] = 0;
  while (X) {
      ++H[0];
      H[H[0]] = X % BASE;
      X /= BASE;
  }
}

Iata, de exemplu, cum se pune pe vector numarul 195:

In sfarsit, initializarea unui vector cu altul se face printr-o simpla copiere (se pot folosi cu succes rutine de lucru cu memoria, cum ar fi FillChar in Pascal sau memmove in C). Pentru eleganta, poate fi folosita si atribuirea cifra cu cifra:

void AtribHuge(Huge H, Huge X) {
  int i;
  for (i = 0; i <= X[0]; ++i) {
    H[i] = X[i];
  }
}

Compararea

Pentru a compara doua numere "uriase", incepem prin a afla numarul de cifre semnificative (deoarece, in urma anumitor operatii pot rezulta zerouri nesemnificative care "atarna" totusi la numarul de cifre). Aceasta se face decrementand numarul de cifre al fiecarui numar atata timp cat cifra cea mai semnificativa este 0. Dupa ce ne-am asigurat asupra acestui punct, comparam numarul de cifre al celor doua cifre. Numarul cu mai multe cifre este cel mai mare. Daca ambele numere au acelasi numar de cifre, pornim de la cifra cea mai semnificativa si comparam cifrele celor doua numere pana la prima diferenta intalnita. In acest moment, numarul a carui cifra este mai mare este el insusi mai mare.
Daca toate cifrele numerelor sunt egale doua cate doua, atunci in mod evident numerele sunt egale.

Dupa cum se vede, algoritmul seamana foarte bine cu ceea ce s-a invatat la matematica prin clasa a II-a (doar ca atunci nu ni s-a spus ca este vorba de un "algoritm"). Rutina de mai jos compara doua numere uriase H1 si H2 si intoarce -1, 0 sau 1, dupa H1 este mai mic, egal sau mai mare decat H2.

int Sgn(Huge H1, Huge H2) {
  // Elimina zero-urile semnificative, daca exista.
  while (H1[0] && !H1[H1[0]]) H1[0]--;
  while (H2[0] && !H2[H2[0]]) H2[0]--;

  if (H1[0] < H2[0]) {
    return -1;
  } else if (H1[0] > H2[0]) {
    return +1;
  }

  for (int i = H1[0]; i > 0; --i) {
    if (H1[i] < H2[i]) {
      return -1;
    } else if (H1[i] > H2[i]) {
      return +1;
    }
  }
  return 0;
}

Adunarea a doi vectori

Fiind dati doi vectori, A cu M cifre si B cu N cifre, adunarea lor se face in mod obisnuit, ca la aritmetica. Pentru a evita testele de depasire, se recomanda sa se completeze mai intai vectorul cu mai putine cifre cu zerouri pana la dimensiunea vectorului mai mare. La sfarsit, vectorul suma va avea dimensiunea vectorului mai mare dintre A si B, sau cu 1 mai mult daca apare transport de la cifra cea mai semnificativa. Procedura de mai jos adauga numarul B la numarul A.

void Add(Huge A, Huge B)
/* A <- A+B */
{ int i,T=0;

  if (B[0]>A[0])
    { for (i=A[0]+1;i<=B[0];) A[i++]=0;
      A[0]=B[0];
    }
    else for (i=B[0]+1;i<=A[0];) B[i++]=0;
  for (i=1;i<=A[0];i++)
    { A[i]+=B[i]+T;
      T=A[i]/10;
      A[i]%=10;
    }
  if (T) A[++A[0]]=T;
}

Scaderea a doi vectori

Se dau doi vectori A si B si se cere sa se calculeze diferenta A - B. Se presupune B ≤ A (acest lucru se poate testa cu functia Sgn). Pentru aceasta, se porneste de la cifra unitatilor si se memoreaza la fiecare pas "imprumutul" care trebuie efectuat de la cifra de ordin imediat superior (imprumutul poate fi doar 0 sau 1). Deoarece B ≤ A, avem garantia ca pentru a scadea cifra cea mai semnificativa a lui B din cifra cea mai semnificativa a lui A nu e nevoie de imprumut.

void Subtract(Huge A, Huge B)
/* A <- A-B */
{ int i, T=0;

  for (i=B[0]+1;i<=A[0];) B[i++]=0;
  for (i=1;i<=A[0];i++)
    A[i]+= (T=(A[i]-=B[i]+T)<0) ? 10 : 0;
    /* Adica A[i]=A[i]-(B[i]+T);
       if (A[i]<0) T=1; else T=0;
       if (T) A[i]+=10; */
  while (!A[A[0]]) A[0]--;
}

Inmultirea si impartirea cu puteri ale lui 10

Aceste functii sunt uneori utile. Ele pot folosi si functiile de inmultire a unui vector cu un scalar, care vor fi prezentate mai jos, dar se pot face si prin deplasarea intregului numar spre stanga sau spre dreapta. De exemplu, inmultirea unui numar cu 100 presupune deplasarea lui cu doua pozitii inspre cifra cea mai semnificativa si adaugarea a doua zerouri la coada. Principalul avantaj al scrierii unor functii separate pentru inmultirea cu 10, 100, ..., este ca se pot folosi rutinele de acces direct al memoriei (FillChar, respectiv memmove). Iata functiile care realizeaza deplasarea vectorilor, atat prin mutarea blocurilor de memorie, cat si prin atribuiri succesive.

void Shl(Huge H, int Count)
/* H <- H*10ACount */
{ 
  /* Shifteaza vectorul cu Count pozitii */
  memmove(&H[Count+1],&H[1],sizeof(int)*H[0]);
  /* Umple primele Count pozitii cu 0 */
  memset(&H[1],0,sizeof(int)*Count);
  /* Incrementeaza numarul de cifre */
  H[0]+=Count;
}

void Shl2(Huge H, int Count)
/* H <- H*10ACount */
{ int i;

  /* Shifteaza vectorul cu Count pozitii */
  for (i=H[0];i;i--) H[i+Count]=H[i];
  /* Umple primele Count pozitii cu 0 */
  for (i=1;i<=Count;) H[i++]=0;
  /* Incrementeaza numarul de cifre */
  H[0]+=Count;
}

void Shr(Huge H, int Count)
/* H <- H/10ACount */
{ 
  /* Shifteaza vectorul cu Count pozitii */
  memmove(&H[1],&H[Count+1],sizeof(int)*(H[0]-Count));
  /* Decrementeaza numarul de cifre */
  H[0]-=Count;
}

void Shr2(Huge H, int Count)
/* H <- H/10ACount */
{ int i;

  /* Shifteaza vectorul cu Count pozitii */
  for (i=Count+1;i<=H[0];i++) H[i-Count]=H[i];
  /* Decrementeaza numarul de cifre */
  H[0]-=Count;
}

Inmultirea unui vector cu un scalar

Si aceasta operatie este o simpla implementare a modului manual de efectuare a calculului. La inmultirea "de mana" a unui numar mare cu unul de o singura cifra, noi parcurgem deinmultitul de la sfarsit la inceput, si pentru fiecare cifra efectuam urmatoarele operatii:

Inmultim cifra respectiva cu inmultitorul;
Adaugam "transportul" de la inmultirea precedenta;
Separam ultima cifra a rezultatului si o trecem la produs;
Celelalte cifre are rezultatului constituie transportul pentru urmatoarea inmultire;
La sfarsitul inmultirii, daca exista transport, acesta are o singura cifra, care se scrie inaintea rezultatului.

Exact acelasi procedeu se poate aplica si daca inmultitorul are mai mult de o cifra. Singura deosebire este ca transportul poate avea mai multe cifre (poate fi mai mare ca 9). Din aceasta cauza, la sfarsitul inmultirii, poate ramane un transport de mai multe cifre, care se vor scrie inaintea rezultatului. Iata de exemplu cum se calculeaza produsul 312 x 87:

Procedura este descrisa mai jos:

void Mult(Huge H, unsigned long X)
/* H <- H*X */
{ int i;
  unsigned long T=0;

  for (i=1;i<=H[0];i++)
    { H[i]=H[i]*X+T;
      T=H[i]/10;
      H[i]=H[i]%10;
    }
  while (T) /* Cat timp exista transport */
    { H[++H[0]]=T%10;
      T/=10;
    }
}

Inmultirea a doi vectori

Daca ambele numere au dimensiuni mari si se reprezinta pe tipul de date Huge, produsul lor se calculeaza inmultind fiecare cifra a deinmultitului cu fiecare cifra a inmultitorului si trecand rezultatul la ordinul de marime (exponentul lui 10) cuvenit. De fapt, acelasi lucru il facem si noi pe hartie. Considerand acelasi exemplu, in care ambele numere sunt "uriase", produsul lor se calculeaza de mana astfel

S-a luat deci fiecare cifra a inmultitorului si s-a efectuat produsul partial corespunzator, corectand la fiecare pas rezultatul prin calculul transportului. Rezultatul pentru fiecare produs partial s-a scris din ce in ce mai in stanga, pentru a se alinia corect ordinele de marime. Acest procedeu este oarecum incomod de implementat. Se pot face insa unele observatii care usureaza mult scrierea codului:

Prin inmultirea cifrei cu ordinul de marime 10i din primul numar cu cifra cu ordinul de marime 10j din al doilea numar se obtine o cifra corespunzatoare ordinului de marime 10i+$j$ in rezultat (sau se obtine un numar cu mai mult de o singura cifra, caz in care transportul merge la cifra corespunzatoare ordinului de marime 10i + j + 1).
Daca numerele au M si respectiv N cifre, atunci produsul lor va avea fie M + N fie M + N - 1 cifre. Intr-adevar, daca numarul A are M cifre, atunci 10^M-1 ≤ A < 10^M si 10^N-1 ≤ B < 10^N, de unde rezulta 10^M+N-2 ≤ A x B < 10^M+N.
La calculul produselor partiale se poate omite calculul transportului, acesta urmand a se face la sfarsit. Cu alte cuvinte, intr-o prima faza putem pur si simplu sa inmultim cifra cu cifra si sa adunam toate produsele de aceeasi putere, obtinand un numar cu "cifre" mai mari ca 9, pe care il aducem la forma normala printr-o singura parcurgere. Sa reluam acelasi exemplu:

Aceasta operatie efectueaza MxN produse de cifre si M+N (sau M+N-1, dupa caz) "transporturi" pentru aflarea rezultatului, deci are complexitatea O(MxN). Iata si implementarea:

void MultHuge(Huge A, Huge B, Huge C)
/* C <- A x B */
{I int i,j,T=0;

  C[0]=A[0]+B[0]-1;
  for (i=1;i<=A[0]+B[0];) C[i++]=0;
  for (i=1;i<=A[0];i++)
    for (j=1;j<=B[0];j++)
      C[i+j-1]+=A[i]*B[j];
  for (i=1;i<=C[0];i++)
    { T=(C[i]+=T)/10;
      C[i]%=10;
    }
  if (T) C[++C[0]]=T;
}

Mai exista o alta modalitate de a inmulti doua numere de cate N cifre fiecare, care are complexitatea $O(N^{log_2 3}) \approx O(N^{1.58}) \approx O(N \sqrt{N})$ . Ea deriva de un algoritm propus de Strassen in 1969 pentru inmultirea matricelor. Diferenta se face simtita, ce-i drept pentru valori mari ale lui N, dar constanta multiplicativa creste destul de mult si, in plus, solutia e mai greu de implementat; de aceea nu recomandam implementarea ei in timpul concursului. Ideea de baza este sa se micsoreze numarul de inmultiri si sa se mareasca numarul de adunari, deoarece adunarea a doi vectori se face in O(N), pe cand inmultirea se face in O(N 2). Sa consideram intregii A si B, fiecare de cate N cifre. Trebuie sa-i inmultim intr-un timp T(N) mai bun decat O(N 2). Sa impartim numarul A in doua "bucati" C si D, fiecare de cate N/2 cifre, iar intregul B in doua bucati E si F, tot de cate N/2 cifre (presupunem ca N este par):

Atunci se poate scrie relatia:
A x B = (C x 10^N/2 + D) x (E x 10^N/2 + F) = CE x 10^N (CF + DE) x 10^N/2 + DF

Pentru a putea calcula produsul A x B avem, prin urmare, nevoie de patru produse partiale, de trei adunari si de doua inmultiri cu puteri ale lui 10. Adunarile si inmultirile cu puteri ale lui 10 se fac in timp liniar. Daca efectuam cele patru produse partiale prin patru inmultiri, rezulta formula recurenta de calcul T(N)=4T(N/2)+O(N) care duce prin eliminarea recurentei la $T(N)\in O(N)^2$ . Cu alte cuvinte, inca n-am castigat nimic. Trebuie sa reusim cumva sa reducem numarul de inmultiri de la 4 la 3, chiar daca prin aceasta vom mari numarul de adunari necesare. Sa definim produsul G = (C + D) x (E + F) = CE + CF + DE + DF = CE + DF + (CF + DE)
Atunci putem scrie: A x B = CE x 10^N + (G - CE - DF) x 10^N/2 + DF
Pentru aceasta varianta, folosim doar trei inmultiri, si chiar daca avem nevoie de sase adunari si scaderi si doua inmultiri cu puteri ale lui 10, complexitatea se va reduce la $O(N^{log_2 3})$ . In cazul in care N este o putere a lui 2, impartirea in doua a numerelor se poate face fara probleme la fiecare pas, pana se ajunge la numere de o singura cifra, care se inmultesc direct. In cazul in care N nu este o putere a lui 2, este comod sa se completeze numerele cu zerouri pana la o putere a lui 2. In functiile descrise mai jos, MultRec nu face decat inmultirea recursiva, pe cand MultVect2 se ocupa si de corectarea numarului de cifre (incrementarea pana la o putere a lui 2). Pentru calculul produselor C x E si D x F, procedura MultRec se autoapeleaza; pentru calcularea produsului (C+D) x (E+F), insa, este nevoie sa fie apelata procedura MultVect2, deoarece prin cele doua adunari poate sa apara o crestere a numarului de cifre al factorilor, care in acest caz trebuie readusi la un numar de cifre egal cu o putere a lui 2.

void MultHuge2(Huge A, Huge B, Huge P);
void MultRec(Huge A, Huge B, Huge P)
{ Huge C,D,E,F,CE,DF;
  if (A[0]==1)
    { P[1]=A[1]*B[1];
      P[0]=(P[2]=P[1]/10)>0 ? 2 : 1;
      P[1]%=10;
    }
    else { P[0]=0;
           AtribHuge(C,A);Shr(C,A[0]/2);
           AtribHuge(D,A);D[0]=A[0]/2;
           AtribHuge(E,B);Shr(E,B[0]/2);
           AtribHuge(F,B);F[0]=B[0]/2;
           MultRec(C,E,CE);MultRec(D,F,DF);
           Add(C,D);Add(E,F);
           MultHuge2(C,E,P);
           Subtract(P,CE);Subtract(P,DF);
           Shl(P,A[0]/2);
           Shl(CE,A[0]);Add(P,CE);
           Add(P,DF);
         }
}
void MultHuge2(Huge A, Huge B, Huge P)
/* P <- A x B, varianta NA(lg 3) */
{ int i,j,NDig=A[0]>B[0] ? A[0] : B[0],Needed=1,SaveA,SaveB;
  /* Corecteaza numarul de cifre */
  while (Needed<NDig) Needed<<=1;
  SaveA=A[0];SaveB=B[0];A[0]=B[0]=Needed;
  for (i=SaveA+1;i<=Needed;) A[i++]=0;
  for (i=SaveB+1;i<=Needed;) B[i++]=0;
  MultRec(A,B,P);
  /* Restaureaza numarul de cifre */
  A[0]=SaveA;B[0]=SaveB;
  while (!P[P[0]] && P[0]>1) P[0]--;
}

Impartirea unui vector la un scalar

Ne propunem sa scriem o functie care sa imparta numarul A de tip Huge la scalarul B, sa retina valoarea catului tot in numarul A si sa intoarca restul (care este o variabila scalara). Sa pornim de la un exemplu particular si sa generalizam apoi procedeul: sa calculam catul si restul impartirii lui 1997 la 7. Cu alte cuvinte, sa gasim acele numere C de tip Huge si $R\in{0, 1, 2, 3, 4, 5, 6}$ cu proprietatea ca 1997 = 7 x C + R.

La fiecare pas se coboara cate o cifra de la deimpartit alaturi de numarul deja existent (care initial este 0), apoi rezultatul se imparte la impartitor ($7$ in cazul nostru). Catul este intotdeauna o cifra si se va depune la sfarsitul catului impartirii, iar restul va fi folosit pentru urmatoarea impartire. Restul care ramane dupa ultima impartire este tocmai R pe care il cautam. Procedeul functioneaza si atunci cand deimpartitul are mai multe cifre. La sfarsit trebuie sa decrementam corespunzator numarul de cifre al catului, prin neglijarea zerourilor create la inceputul numarului. Numarul maxim de cifre al catului este egal cu cel al deimpartitului.

unsigned long Divide(Huge A, unsigned long X)
/* A <- A/X si intoarce A%X */
{ int i;
  unsigned long R=0;

  for (i=A[0];i;i--)
    { A[i]=(R=10*R+A[i])/X;
      R%=X;
    }
  while (!A[A[0]] && A[0]>1) A[0]--;
  return R;
}

Daca dorim numai sa aflam restul impartirii, nu mai avem nevoie decat sa recalculam restul la fiecare pas, fara a mai modifica vectorul A:

unsigned long Mod(Huge A, unsigned long X)
/* Intoarce A%X */
{ int i;
  unsigned long R=0;

  for (i=A[0];i;i--)
    R=(10*R+A[i])%X;
  return R;
}

infoarena informatica de performanta