Introduction à MPIR

Thibaut Cuvelier

MPIR, pour Multiple Precision Integers and Rationals, est, comme son nom l'indique, une bibliothèque de calcul sur des entiers et des rationnels à précision multiple (arbitraire). Elle est développée, à en croire le site officiel, par des mathématiciens et des informaticiens, tant amateurs que professionnels.

Il s'agit d'un fork de GMP (GNU Multiple Precision), ayant comme but d'avoir une communauté developer- and vendor-friendly (une des raisons du fork a été le refus actif de l'équipe de GMP de supporter les plateformes Microsoft).

Commentez

- Titre : Introduction à MPIR

- Auteur : Thibaut Cuvelier

- Parution : 29 juillet 2011

- Mise à jour : 7 août 2011

- Licence : Copyright ® 2011 Thibaut Cuvelier. Aucune reproduction, même partielle, ne peut être faite de ce site et de l'ensemble de son contenu : textes, documents, images, etc. sans l'autorisation expresse de l'auteur. Sinon vous encourez selon la loi jusqu'à trois ans de prison et jusqu'à 300 000 € de dommages et intérêts.

Sommaire

- I. Deux mots sur MPIR

- II. L'installation

- III. Les bases de MPIR

- III-A. Les en-têtes

- III-B. Les types

- III-C. Les conventions

- III-D. L'initialisation des variables

- IV. Les entiers

- IV-A. Une variable

- IV-B. Son initialisation

- IV-C. Son affectation

- IV-D. Son initialisation et son affectation en même temps

- IV-E. Sa conversion

- IV-F. Des opérations arithmétiques

- V. Les rationnels

- VI. Les nombres en virgule flottante

- VII. MPIR et le C++

- VIII. Conclusion

I - Deux mots sur MPIR

MPIR est une bibliothèque portable écrite en C pour une arithmétique à précision arbitraire sur des entiers, des rationnels et des nombres à virgule flottante, voulant être la plus rapide pour ceux qui sont limités par les types de base du C. Une API C++ est elle aussi disponible.

Cette bibliothèque est prévue autant pour ceux qui traitent de petits nombres (quelques centaines de bits) que pour ceux qui ont des données de plus grande ampleur (des milliers, voire des millions de bits par nombre) : les algorithmes s'adaptent en fonction de la taille des entrées pour donner de bonnes performances dans tous les cas. (1) La rapidité est largement préférée à la simplicité et à l'élégance du code.

Cette introduction ne veut pas être une introduction en profondeur aux concepts sous-jacents à MPIR, à l'arithmétique sur bignums, simplement une introduction à l'utilisation de base de MPIR.

On utilisera ici MPIR en version 2.4.0 ; cependant, étant donné que l'on ne s'avance pas loin dans son API, il devrait rester exploitable avec les versions antérieures et postérieures de MPIR, ainsi qu'avec GMP.

II - L'installation

Tout d'abord, il faut vérifier que des binaires ne sont pas disponibles pour votre plateforme ; certaines distributions Linux en disposent, mais ce ne sont pas toujours les versions les plus récentes. Vérifiez dans votre gestionnaire de paquet paquet de quoi il retourne.

Sinon, il faudra télécharger les sources de MPIR et les compiler. Pour ce faire, direction la page d'accueil de la bibliothèque, section Downloads, téléchargez la dernière version (la version LZMA est plus légère car mieux compressée que la BZ2). Il faut ensuite décompresser ces sources dans un endroit de travail (peu importe, tant qu'aucun fichier existant n'est écrasé dans la partie, pour éviter tout problème ultérieur).

Ensuite, avec les outils de compilation GNU (généralement avec MinGW ou GCC), la manipulation n'est pas compliquée :

./configure --prefix=/usr/local --enable-fat --enable-cxx --enable-gmpcompat
make
make install

Le paramètre --prefix permet de choisir l'endroit d'installation (la valeur présentée ci-dessous est celle par défaut).

--enable-fat, quant à lui, demande de compiler des binaires plus lourds que la moyenne (uniquement pour les architectures x86 et x86_64, les plus répandues), car ils contiennent plusieurs versions du code assembleur, la plus rapide étant sélectionnée à l'exécution.

Finalement, --enable-cxx demande de compiler le binding C++, un compilateur C++ est cependant requis. Cela générera un nouveau fichier d'import des bibliothèques spécifique au C++. De même, --enable-gmpcompat demandera la compilation de la couche de compatibilité avec GMP, ce qui permettra d'utiliser MPIR en parfait remplaçant de GMP, sans la moindre modification : ainsi, de nombreux projets tels que MPFR, ECM, MPC ou PARI pourront utiliser MPIR sans subir de modification manuelle.

Beaucoup d'autres options sont disponibles, n'hésitez pas à regarder le contenu de l'aide à ce sujet.

Pour les utilisateurs de Visual Studio 2010, des fichiers de projet sont disponibles dans le dossier build.vc10 de la distribution de MPIR ; consultez le fichier d'aide README dans ce dossier, la compilation pouvant être configurée autant que précédemment, mais d'une manière différente.

III - Les bases de MPIR

III-A - Les en-têtes

En C, sans en-tête, pas de répit. On commence donc par inclure l'en-tête C de MPIR. Si des fonctionnalités basées sur des fichiers sont requises, il faut inclure stdio.h avant l'en-tête de MPIR ; pour des prototypes avec paramètres va_list, stdarg.h ou varargs.h ; pour des fonctions avec des paramètres struct obstack, obstack.h. Pour bénéficier de toutes les fonctionnalités, on fera donc :

#include <stdio.h>
#include <stdarg.h>
#include <varargs.h>
#include <mpir.h>

Cependant, ces fonctionnalités ne seront pas abordées dans ce tutoriel, cette simple inclusion suffira largement :

#include <mpir.h>

Pour que la compilation se passe ensuite sans souci, il faut encore préciser les bibliothèques à lier (dynamiquement ou statiquement) en ligne de commande. Pour GCC, cela donne :

gcc myprogram.c -lmpir

Le même genre d'options sont disponibles pour les autres suites d'outils de compilation, leur documentation reste une référence à ce sujet.

III-B - Les types

MPIR est concentré autour de trois types de base :

- mpz_t : un entier signé ;

- mpq_t : un rationnel (quotient de deux entiers) ; (2)

- mpf_t : un nombre en virgule flottante (mantisse de précision arbitraire mais exposant de précision limitée).

À chaque fois, on ne garde qu'un pointeur sur les objets avec les champs requis par MPIR. Ces objets sont au final assez légers : les tailles et les pointeurs sur les espaces alloués pour les données. L'allocation de mémoire est gérée automatiquement : dès que plus d'espace est requis, il est alloué automatiquement.

mpz_t et mpq_t ne peuvent que grandir en mémoire pour éviter de fréquentes réallocations (3) ; a contrario, mpf_t utilise une quantité d'espace prédéfinie par la précision désirée, la taille d'une telle variable ne variera donc pas.

III-C - Les conventions

On peut remarquer quelques points communs : tous les types correspondent au masque mp?_t. t signifie qu'il s'agit d'un type ; le préfixe mp? indique le type de nombres géré. Ainsi, on retrouvera ce même préfixe pour les fonctions qui gèrent ces types : mpz_add, par exemple, additionne deux entiers.

Revenons sur cette fonction d'addition : quels sont ses paramètres ? Elle prend trois opérandes de type mpz_t : le résultat et les deux opérandes de l'addition. Dans MPIR, les fonctions ont toujours cet ordre de paramètres : d'abord les sorties, puis les entrées, par analogie avec l'opérateur d'assignation. Les opérandes d'entrée sont généralement définies const : elles ne seront en aucun cas modifiées. Il est donc facilement possible pour une fonction de retourner plusieurs variables. On remarquera que cela laisse le retour de la fonction (par l'instruction return) pour, par exemple, vérifier qu'une opération s'est bien déroulée.

III-D - L'initialisation des variables

Toute variable d'un type déjà cité doit être initialisée avant toute utilisation ; on utilisera la famille de fonctions mp?_init pour ce faire. Une fois initialisée, une variable peut recevoir une valeur autant de fois que nécessaire (être assignée).

Aussi, il n'est pas possible de donner une valeur à une variable de type MPIR à l'aide du simple opérateur =, cela n'étant possible que pour les types de base (la surcharge des opérateurs n'existant pas en C).

Le couple d'opérations initialisation-assignation étant très fréquent, des fonctions permettent de combiner les deux opérations : cela permet de gagner un peu de temps à l'exécution (chaque appel de fonction étant une perte de temps potentielle ; de plus, il est possible d'optimiser encore plus ces opérations).

IV - Les entiers

IV-A - Une variable

On commence toujours par définir une variable :

mpz_t i;

IV-B - Son initialisation

Ensuite, on doit l'initialiser (la valeur par défaut est 0) :

mpz_init(i);

On peut, si l'on a une estimation de la taille requise pour cette variable, la définir à l'initialisation, ce qui aura pour effet d'éviter de coûteuses réallocations :

mpz_init2(i, 15684);

On ne définit ici qu'une taille initiale de la variable, elle continuera à grandir avec les besoins.

La gestion mémoire a toujours été importance en C, il est donc possible (et plus que chaudement recommandé) de libérer la mémoire allouée :

mpz_clear(i);

IV-C - Son affectation

Diverses méthodes d'affectation sont possibles : on peut copier un entier déjà affecté ou bien utiliser un nombre dans un type standard du C, on peut également utiliser une chaîne (un tableau de char, la fin de la chaîne est indiquée par \0), auquel cas on peut spécifier la base de l'entrée (de deux à trente-six sans distinction entre majuscules et minuscules, les lettres étant utilisées comme chiffres ; entre trente-sept et soixante-deux, les majuscules représentent les chiffres entre dix et trente-cinq tandis que les minuscules terminent la ronde entre trente-six et soixante-et-un).

Fonction

Type du second opérande

mpz_set

mpz_t

mpz_set_ui

unsigned long int

mpz_set_si

signed long int

mpz_set_ux

uintmax_t

mpz_set_sx

intmax_t

mpz_set_d

double

mpz_set_q

mpq_t

mpz_set_f

mpf_t

Ces fonctions ont un prototype semblable à void mpz_set*(mpz_t rop, type op) : une sortie et une entrée, le contenu de op est mis dans rop.

Pour des chaînes, on peut utiliser int mpz_set_str(mpz_t rop, char * str, int base) (les espaces blancs sont ignorés), qui retourne 0 si la chaîne est un nombre valide dans la base donnée, -1 sinon.

On peut aussi échanger les valeurs de deux variables de manière efficace avec void mpz_swap(mpz_t rop1, mpz_t rop2).

IV-D - Son initialisation et son affectation en même temps

On ajoute simplement le mot-clé init dans les fonctions d'affectation ; il est inutile de définir la taille voulue, étant donné qu'elle est déterminée avec l'opérande d'entrée. Exemple complet :

#include <mpir.h>

int main(void)
{
mpz_t i;

mpz_init_set_ui(i, (unsigned long int) 42);

mpq_clear(i);

return 0;
}

IV-E - Sa conversion

Pour transmettre un résultat, il n'est pas toujours aisé d'utiliser un type de MPIR, il est donc possible de revenir aux types standard du C : exception, un seuleopérande (la variable à convertir), la sortie se fait par return (plus pratique si l'on doit, par exemple, afficher le résultat par printf()).

La bibliothèque est cohérente : on peut reprendre le tableau précédent pour les types, il suffit de changer le radical utilisé par mpz_get.

Il faut cependant faire attention lors de la conversion : MPIR est parfois amené à effectuer des opérations sur des nombres extrêmement grands, potentiellement impossibles à représenter dans les types standard du C. En fonction des cas, le signe est ignoré (type de destination unsigned) ou bien seuls les bits de poids faible (autant que possible) sont renvoyés.

Heureusement, on dispose de la fonction mpz_fits_slong_p() : elle permet de déterminer si la variable peut être convertie sans perte en signed long int.

Il faut aussi noter la fonction mpz_sizeinbase (op, base), qui permet de déterminer la taille requise pour afficher le nombre dans la base donnée. Avec deux octets supplémentaires (pour le signe et \0), cela donne la taille maximale d'une chaîne de caractères représentant le nombre ; on obtient cette chaîne de caractères avec la fonction char * mpz_get_str (char * str, int base, mpz_t op). Le paramètre char * str n'est pas très naturel : s'il est NULL, alors MPIR allouera la mémoire nécessaire pour afficher le nombre ; sinon, il utilisera ce pointeur pour commencer à stocker la représentation du nombre.

IV-F - Des opérations arithmétiques

La partie la plus intéressante : on peut effectuer des opérations arithmétiques avec MPIR ! Vu leur extrême complexité, elles sont présentées dans un tableau, avec la correspondance en termes d'opérateurs C. D'autres versions sont disponibles et documentées dans le manuel (on peut, par exemple, remplacer certains opérandes par des types de base du C, sans besoin d'utiliser les types de MPIR).

Fonction

Équivalent en (pseudo-)C

void mpz_add(mpz_t rop, mpz_t op1, mpz_t op2)

rop = op1 + op2;

void mpz_sub(mpz_t rop, mpz_t op1, mpz_t op2)

rop = op1 - op2;

void mpz_mul(mpz_t rop, mpz_t op1, mpz_t op2)

rop = op1 * op2;

void mpz_addmul(mpz_t rop, mpz_t op1, mpz_t op2)

rop = rop + (op1 * op2);

void mpz_submul(mpz_t rop, mpz_t op1, mpz_t op2)

rop = rop - (op1 * op2);

void mpz_neg(mpz_t rop, mpz_t op)

rop = - op;

void mpz_abs(mpz_t rop, mpz_t op)

rop = fabs(op);

void mpz_powm(mpz_t rop, mpz_t base, mpz_t exp, mpz_t mod)

rop = exp(base, exp) modulo mod;

Les fonctions de division sont, quant à elles, beaucoup plus variées : on peut toujours arrondir vers le haut, le bas ou tronquer le résultat, notamment. Elles sont suffisamment abondantes pour qu'il soit inutile de les lister ici. De même pour l'extraction des racines, la détection des nombres premiers probables, le PGCD, le PPCM, la factorielle et bien d'autres. De plus, ces fonctions n'intéresseront pas forcément tous les publics tout en étant correctement documentées. Pour information : il existe aussi des fonctions pour comparer, inverser, calculer les symboles de Jacobi ou de Legendre, calculer un terme de la suite de Fibonacci, effectuer des opérations logiques et binaires, lire ou écrire le résultat dans un fichier, générer des nombres aléatoires.

V - Les rationnels

On récupère les mêmes fonctions de base que précédemment :

#include <mpir.h>

int main(void)
{
mpq_t i;
mpq_init(i);
mpq_set_ui(i, (unsigned long int) 42, (unsigned long int) 7); // i vaut maintenant 42/7

mpq_t j;
mpq_init(j);
mpq_set_ui(j, (unsigned long int) 87, (unsigned long int) 7);
mpq_canonicalize(j); // On s'assure que j est bien simplifiée

mpq_sub(i, i, j); // i == i - j   ou   i -= j

mpz_add(mpq_numref(i), mpq_denref(i), mpq_numref(j)); // Plus intéressant : on peut accéder
// séparément au numérateur et au dénominateur. mpq_numref() renvoie le mpz_t du numérateur
// de son argument ; mpq_denref(), le dénominateur.

mpq_clear(i);
mpq_clear(j);

return 0;
}

VI - Les nombres en virgule flottante

Ces nombres sont plus particuliers en ce que l'on peut utiliser la précision par défaut ou celle précisée, leur taille ne changera jamais. Sinon, ils sont manipulables comme précédemment, en utilisant des fonctions appropriées.

#include <mpir.h>

int main(void)
{
mpf_t x, y;
mpf_init (x); // Précision par défaut
mpf_init2 (y, 2048); // Précision d'au moins 2048 bits ; on peut récupérer la précision avec
// mpf_get_prec(x), la définir à au moins p bits avec mpf_set_prec(x, p) avec troncature
// ou sans troncature avec mpf_set_prec_raw(x, p).

mpf_set_si(x, (signed long int) (- 42 / 3));
mpf_set_d(y, (double) M_PI);

pdf_add(x, x, y) ; // x += y

printf("%f", mpf_get_d(x));

mpf_clear (x);
mpf_clear (y);

return 0;
}

Bien que disponible, ce module doit être considéré comme obsolète car très peu testé. Les développeurs de MPIR conseillent d'utiliser MPFR à la place.

VII - MPIR et le C++

Un binding C++ est disponible avec MPIR de base si l'on en demande la compilation. Il faudra inclure le fichier mpirxx.h pour en disposer (inutile d'inclure mpir.h), ainsi que de lier le programme aussi avec mpirxx :

g++ mycxxprog.cc -lmpirxx -lmpir

Trois classes sont définies : mpz_class, mpq_class et mpf_class. Grâce aux constructeurs, il n'est pas nécessaire d'initialiser ces variables, cela est fait automatiquement ; avec l'aide de la surcharge des opérateurs, on peut effectuer une affectation avec le signe d'égalité = et utiliser les opérateurs arithmétiques du langage (+, -, * notamment) ; on peut aussi utiliser ces objets directement dans des flux.

#include <iostream>
#include <mpirxx.h>

int main(voic)
{
mpz_class x;
mpz_class y;
mpz_class z;

x = 42; // depuis un int
y = "-84"; // depuis une chaîne de caractères
z = x + 2 * y; // cette expression utilise des variables temporaires ; par contre, z = x + y n'en a pas besoin

std::cout << x << " + 2 * " << y << " = " << z << std::endl;

return 0;
}

On peut récupérer les types C sous-jacents avec les méthodes get_mpz_t, get_mpq_t et get_mpf_t.

VIII - Conclusion

MPIR n'est pas une bibliothèque très compliquée à utiliser, loin de là, et permet d'obtenir des performances quand cela est réellement requis par l'application. Il faut également remarquer que rien n'est laissé au hasard dans les algorithmes, ils sont généralement éprouvés, mais surtout documentés pour ce projet, ce qui facilitera la compréhension du code interne de la bibliothèque. On remarque aussi, en lisant cette section, que des macros sont définies pour savoir quelle fonction utiliser sur la plateforme en cours en fonction des tailles des entrées (par exemple, entre une fonction pour calculer directement un carré ou la multiplication par soi-même, cela pouvant arriver si la fonction est disponible nativement sur le processeur).

MPIR n'est pas à délaisser au profit de nouvelles technologies comme le GPGPU : ses développeurs comptent évoluer dans le sens d'une plus grande parallélisation du code, notamment à l'aide de CUDA et d'OpenMP.

Merci à Claude Leloup pour sa relecture orthographique ! Merci à Jason de l'équipe MPIR pour sa relecture technique !

(1) Au contraire de GMP, qui utilise des algorithmes plutôt asymptotiques : plus l'entrée est de grande taille, plus efficace sera l'algorithme déployé.

(2) Les rationnels sont automatiquement simplifiés (pas de facteur commun entre le numérateur et le dénominateur), car cela minimise la taille des opérandes pour les calculs à venir et donc maximise la vitesse.

(3) Pour éviter trop de réallocations en cours de route, il est possible de définir l'espace requis (en bits) pour la variable si une telle estimation est possible.