Mathématiques et grands nombres

 

Tous ceux qui, un jour, ont eu à rendre un programme, ou pire à faire pour leur propre satisfaction, un programme utilisant des fonctions mathématiques utilisant de grands nombres, se sont peut-être heurtés au problème des limitations de type du Turbo Pascal. Grande nouveauté avec le BP7, les formats COMP, mais les nombres sont encore très limités. Suite à un projet de recherche de grands nombres premiers, qui sont tous positifs, car je ne parlerai que de grands nombres positifs, la première chose que mon binôme à faite, a été d'utiliser les listes, dont il avait l'habitude. Mais tout ceux qui sont habitués à utiliser le Pascal sur des stations VAX/VMS se retrouvent marris devant leurs pauvres PC quand ils s'aperçoivent de la terrible limitation de leur système en ressource dynamique : le Pascal ne désallouant pas après un appel (j'envoie 2 listes, et j'en renvoie une dans la 1e i.e. NA := add(NA,NB), l'ancien NA reste en mémoire.

Dans ce cas-là, que faire ? Passer en mode protégé est une solution, mais quand on utilise des routines asm, tout recalculer c'est la poisse..., ou tout bêtement créer de nouvelles fonctions de grands nombre. C'est donc cette méthode que j'ai choisie pour contourner cet écueil, et dont la seule limitation était la mémoire dynamique du système (on peut toujours y appliquer le mode protégé). L'idée est simple, deux tableaux en l'occurrence pour une suite de fibonnacci, tous deux de 32000 cases de byte, puis créer de nouvelles fonctions d'addition, de soustraction ... Avantage immédiat : une case de liste fait 5 byte, une de tableau 1, d'où un gain conséquent de mémoire.

Le type est défini juste ici :

        type
          numtab = array[1..32000] of byte;

        BIGnum = record
          t : ^numtab;
          s : integer;
        end;

Le record Bignum a deux champs, le premier T le pointeur sur tableau de byte, et le deuxième, la taille du nombre contenu dans T, très utile pour optimiser le calcul.

Pour les utiliser, on les crée par :

        var NA : BIGnum; (* on les définit *)

        begin
          new(NA.T); (crée ce nombre, A NE JAMAIS OUBLIER SOUS PEINE DE MORT)
          init_BN(NA);   (on le met a 0)

          add_c(NA,2);   (on met 2 dans NA, c.f. plus loin...)

          release(NA.T); (détruit NA)
        end;

Voyons les fonctions mathématiques de base :

    {------[ Procédures mathématiques ]-------------------------17/05/95-}
On initialise en mettant des zéro partout et en mettant la taille à 0 :

        procedure init_BN(var ptr:BIGnum);
                                        (Initialise un grand nombre à zéro)
        var i : integer;
        begin
          for i := 1 to 32000 do
            ptr.t^ [i] := 0;
          ptr.s := 0;
        end;

Renvoyer le plus grand est souvent nécessaire pour deux Integers :

        function max(i,j:integer):integer;
                                       (Renvoie le plus grand des Integers)
        begin
          if i>j then
            max := i
          else
            max := j;
        end;

Il faut déjà implémenter une fonction pour savoir si on a affaire à un nombre nul, on ne sait jamais, mais ça va servir, je vous le promets :

        function est_nul(NA:BIGnum):boolean;
        begin
          est_nul := (NA.s=0) OR ( (NA.s=1) AND (NA.t^[1]=0) )
        end;

Première procédure intéressante après l'initialisation, l'addition de base, c'est à dire celle de plus bas niveau, l'addition par une constante inférieure à 10. (C.f. cours de numérologie ). J'appelle par la suite BN un Big Number, un grand nombre.

        procedure add_c(var NA:BIGnum;cst:byte);
        var ret : byte;           (la retenue)
            i   : integer;        (le compteur)
            c   : integer;        (le calcul intermédiaire)
        Begin

          ret := 0;                 (pas de retenue au départ)

        if not((cst=0) or est_nul(NA)) then         (zéros, rien a traiter)
          for i := 1 to na.s+1 do     (on va de 1 a n+1, en cas de retenue)
          begin

            if i>1 then cst := 0;        (idiot, mais utile pour optimiser)

            c := NA.t^[i] + (cst + ret);      (valeur actuelle plus cst si)
                                              (premier rang, et retenue si)
                                              (il y a lieu)
            ret := 0;                       (on l'a utilisée, on la retire)

            if (c>=10) then           (si il y a retenue dans notre calcul)
            begin
              ret := 1;                   (on la met pour le prochain tour)
              c := c-10;                  (on retire 10 de c)
            end;

            NA.t^[i] := c;     (on met dans la position du BN, la valeur c)

            if (ret=0) then               (plus de retenue, on peut sortir)
              exit;                       (allez, hop on sort, c'est fini)
                                        
          end;                            (fin du FOR BEGIN)

          if (NA.t^[na.s+1]<>0) then inc(NA.s);         (si le BN a grandi)

        end;

Et maintenant l'addition complète :

        procedure add(var NA, NB, ND : BIGnum);
        var ret : byte;
            i : integer;
            c : integer;
        Begin
          ret := 0;                 (pas de retenue)

          if not(est_nul(NB) or est_nul(NA)) then          (dès 0, on sort)
          begin                                   
            for i := 1 to max(NA.s,NB.s)+1 do    (on calcule au plus grand)
                                                 (des deux+1 si retenue)
            begin
              c := NA.t^[i] + NB.t^[i] + ret;
              ret := 0;
              if c>=10 then
              begin
                ret := 1;
                c := c - 10;
              end;
              ND.t^[i] := c;
            end;
            if (ND.t^[i]<>0) then ND.s := i     (si il y a du monde en n+1)
            else if (i>1)and(ND.t^[i-1]<>0) then ND.s := i - 1;     (sinon)
          end
          else
            ND:=NA;
        end;

Rien de bien nouveau et pas d'astuces, mais voyons la soustraction d'une constante :

        procedure sous_c(var NA : BIGnum;cst : byte);
        var ret : byte;
            i : integer;
            c : integer;
        Begin

          ret := 0;                                        (pas de retenue)

          if not est_nul(NA)and(cst<>0) then  (si NA est nul, rien a faire)
            for i := 1 to NA.s do

            begin

              if i>1 then cst := 0;                            (c.f. add_c)

              if (NA.t^[i]<(cst+ret)) and (i0) and (NA.t^[i]=0) then               
    		    dec(NA.s);                       (si pas dernière case)

                    exit;

              end
            end;

            if (NA.t^[i]<>0) then NA.s := i           (si taille inchangée)
            else if (i>1) and (NA.t^[i-1]<>0) then NA.s := i-1;     (sinon)

        end;

Ici une optimisation, renvoie le plus grand, comme on travaille seulement en nombres positifs, donc on soustrait le plus grand au plus petit :

       
 function greater(NA,NB:BIGnum):boolean;
        var ret : byte;
            i : integer;
            c : integer;
        Begin

          if NA.s>NB.s then            (on commence par regarder la taille)
            greater:=true
          else 
	    if NA.s<NB.s then
              greater:=false
            else                         (et après comparaison case a case)
            begin
              i:=max(na.s,nb.s);                               (taille max)
              while (NA.t^[i]=NB.t^[i])and(i>0) do   (boucle tant qu'égaux)
                dec(i);               (sens décroissant, c'est plus rapide)
              if (NA.t^[i]>=NB.t^[i]) then
                greater:=true
              else
                greater:=false;
            end;

        end;

Et maintenant la soustraction, sous vos applaudissements :

        procedure sous(var NA,NB,ND:BIGnum);
        var ret : byte;
            i : integer;
            c : integer;
        Begin
          ret := 0;
          if est_nul(NB) then
            ND:=NA
          else
            for i := 1 to max(na.s,nb.s) do                   (le parcours)
            begin

              if (NA.t^[i]< (NB.t^[i]+ret) ) then             (A<B+retenue)
              begin

                ND.t^[i] := 10 + NA.t^[i] - ( NB.t^[i] + ret );
                ret := 1;

              end
              else                                            (A>B+retenue)
              begin

                ND.t^[i] := NA.t^[i] - ( NB.t^[i] + ret );
                ret := 0;

              end;

            end;

            if (ND.t^[i]<>0) then ND.s := i
            else if (i>1) and (ND.t^[i-1]<>0) then ND.s := i - 1;
        end;

Voyons maintenant d'autres fonctions plus diverses, par exemple pour savoir si il y a reste dans la division de A par B :

        function modulozero(var NA,NB:BIGnum):boolean;
        var naa : BIGnum;
        Begin
          new(NAa.t);                                (on crée une nouvelle)
          init_BN(NAa);                                   (on l'initialise)
          add(NA,NAa,NAa);                                  (on met NA=NAa)

          while greater(NAa,NB) do    (puis on soustrait NB a Naa tant que)
            sous(NAa,NB,NAa);         (Naa>NB)

          modulozero:= not (is_void(NAa));

          dispose(NAa.t);                     (on nettoie avant de quitter)
        end;

Et maintenant la division par 2, mon cours de numérologie me dit que diviser un nombre par 2, c'est trouver le nombre de 2 dans celui-ci, je vais donc soustraire A div 2 fois 2 pour trouver A div 2 :

        procedure div2(var NAA:BIGnum);
        var t : BIGnum;
            m : longint;
            i : integer;
        Begin

          m := 0;
          if (naa.s<=10) and (naa.t^[10]<=2) then
          begin      (taille < 10, on peut utiliser des longint, compat. tp6)
                 
            for i := naa.s downto 1 do
              m := m*10 + naa.t^[i];               (on passe de BN a longint)
            m := m div 2;                                   (on divise par 2)
                                          
            init_BN(naa);                          (on passe de longint a BN)
            i := 0;
            while (m<>0) do
            begin
              inc(i);
              naa.t^[i] := m mod 10;                    (on prend le dernier)
              m := m div 10;                                (et on le retire)
            end;
            naa.s := i;                                 (on ajuste la taille)
          end
          else                  (on est condamné à des calculs très longs...)
          begin                
            new(t.t);
            init_BN(t);

            repeat
              sous_c(NAA,2);                                    (on retire 2)
              add_c(t,1);                             (et on le comptabilise)
            until (NAA.s<=2) and (NAA.t^[1]<=1) and (NAA.t^[2]=0);

            init_BN(NAA);
            add(NAA,t,NAA);                  (on renvoie t qui est naa div 2)

            dispose(t.t);
          end;

        end;

Et après ca, il faut les afficher ces fichus nombres :

        procedure show(ptr:BIGnum);
        var i : integer;
            s : string;
        begin

          for i := ptr.s downto 1 do
          begin
            str(ptr.t^[i],s);
            affiche(s);        (on affiche un digit, une fonction de votre cru)
          end;

        end;

 

Limitations :

Evidemment, on ne travaille qu'avec des nombres positifs, mais dans le record BIGnum, il est très facile de rajouter un booléen de signe. La mémoire, chaque tableau de 32000 cases fait ... 32 Ko au moins, penser à allouer une Heap conséquente. {$M , , } ...

Il y a des bugs, ce n'est pas optimisé... je vous vois venir, le but de cet article est de présenter des méthodes de grands nombres, si vous avez mieux, ca m'intéresse, et j'attends vos critiques.

 

Conclusion :

Vous voici tout clef en main, vous n'avez plus qu'à faire du copier-coller ... L'exemple est le programme FIB.EXE.

Bibliographie :

Dr Dobb's #226 janvier 95, avec un article sur les grands nombres en C++.

Principes fondamentaux de l'informatique, Ullman ed.DUNOD que je vous recommande vivement.

Remerciements :

Mon prof A. Souah pour ses supers cours, bien que j'attends toujours la dernière éval et mon cours de numération,

HP Mâd, pour son cours d'assembleur Flash qui me « sert » de modèle,

Murdoc, mon binôme, qui malgré ses défauts n'en est pas moins sympa,

et Ecureuil pour ses conseils, et ses routines sans qui certains de mes programmes ne seraient pas ce qu'ils sont.

Et toi, hypocrite lecteur, mon ami, mon frère. ( Baudelaire )