Pari/GP

Posted on by
5

Alors que je me promenais sur la page consacrée à sa classe de MPSI, je suis tombé sur un petit document de PB présentant brièvement le langage Python pour faire des calculs sur ordinateur… Ça m’a rappelé (avec émotion je dois avouer) mes débuts en programmation, lorsque je faisais des programmes php pour calculer la suite de Fibonnaci, factoriser un nombre, etc.

Je voulais donc profiter de ce message pour vous présenter un fabuleux système de calcul formel, PARI/GP, qui est principalement utilisé en théorie des nombres mais qui contient beaucoup de fonctions utiles en mathématiques. De plus, il est libre et distribué selon les termes de la licence publique générale GNU. Il est développé par des français de l’université de Bordeaux I.

PARI/GP est avant tout destiné aux théoriciens des nombres, mais peut-être utilisé avec profit dans d’autres cas (algèbre linéaire, arithmétique, etc.). Il est composé de deux parties :

  • PARI : bibliothèque écrite en C, mais qui peut-être appelée d’un langage haut niveau (Python, C++, Perl, etc.)
  • gp : interface en ligne de commande permettant l’accès aux fonctions de PARI. Le langage utilisé ressemble à du C.

Les calculs sont effectués en précision arbitraire. Il est également très utilisé en cryptologie (vérification qu’un nombre est premier via l’algorithme AKS, calculs sur les courbes elliptiques, algorithme LLL, etc.). En revanche, ses capacités d’intégration et de dérivation sont moins efficaces que des logiciels de calcul formel comme Xcas, Maple ou Mathematica. Il a également une fonction « graphique » qui me semble limitée (mais j’avoue ne pas avoir regardé si il y avait une possibilité d’interaction avec gnuplot ou autre.).

Mieux que de faire encore de nombreux discours, voici quelques exemples d’utilisation. Dans votre terminal, tapez gp pour accéder à l’interface GP pour PARI.

« Hello World » en PARI/GP

Dans tous les langages de programmation, on commence par une telle application, donc je ne déroge pas à la règle:
print("Hello World")
Aucune difficulté, aucun commentaire.

Une boucle en PARI/GP

Montrons la syntaxe d’une boucle for :
for(i=0, 5, print(i))
Ceci affiche les entiers de 0 à 5 à la ligne.

Pour avoir de l’aide sur une fonction, on peut taper ?for et on obtient :
for(X=a,b,seq): the sequence is evaluated, X going from a up to b.

De même, lorsque l’on cherche une fonction dont on ne connaît pas exactement le nom, on peut chercher avec trois points d’interrogations. Exemple :

???matker

retourne

bnfsignunit idealintersect matinverseimage matker
matkerint

Calcul avec des nombres premiers

La fonction isprime est une fonction déterministe (via l’algorithme AKS) qui détermine si un nombre est premier. Jusquà 500 bits de longueur, le calcul est encore rapide pour la primalité :
isprime(2^(501)-1)
retourne 1 en 310 ms environ (pour obtenir le temps de calcul, taper ## après votre dernière opération), alors qu’il met plus de 5 secondes pour un nombre premier de 1024 bits.

D’ailleurs une manière facile d’obtenir un nombre premier aléatoire est d’utiliser la fonction nextprime qui donne le nombre premier qui suit l’argument donné. Par exemple

nextprime(2^(1024)+random(2^(1024)-1))

renvoit le nombre

179769313486231590782518276703737196373630355166285016700028966421710165488969543617
6526889863117584648152455111700332203590123411794386808044158416082415271352025604694004045
5111003982745778025452513076910528815189617790566856810276323839016028883085338806026638240
6524494940761227538110799658431716130780389

qui est un nombre premier de 1024 bits.

Calculs matriciels

On peut faire facilement avec PARI/GP de l’algèbre linéaire. Par exemple :
C=matid(5); C[1,3]=9; C[2,5]=6; b=[1,2,0,1,1];
ne renvoit rien (à cause des ;) mais a bien enregistré la matrice C que l’on peut afficher :
[1 0 5 0 0]
[0 1 0 0 6]
[0 0 1 0 0]
[0 0 0 1 0]
[0 0 0 0 1]

et le vecteur b.

On peut alors demander de résoudre l’équation Cx=b d’inconnue x :

a=matsolve(C, b~)

qui renvoit

[1, -4, 0, 1, 1]~

Et on peut vérifier par

C*a-b~

qui renvoit bien le vecteur nul.

Notez que le signe ~ représente la transposition, et que PARI/GP différentie les vecteurs colonnes et lignes. On peut également accéder à la première ligne de la matrice C[1,] ou sa cinquième colonne C[,5].

Procédures en PARI/GP

On peut bien évidemment faire plusieurs opérations en même temps avec des procédures. Par exemple, on peut mettre le code suivant dans un fichier chiffres.gp (dans le dossier où vous êtes avec le terminal, ou on lui donne l’adresse complète)
Chiffres(d,l,n)=
{
local(i,m,v);

if(type(d) != "t_INT" || d <= 0, print("d doit être un entier strictement positif"); return([]));
if(type(l) != "t_INT" || l <= 0, print("l doit être un entier strictement positif"); return([]));
if(type(n) != "t_INT" || n < 0, print("n doit être un entier positif"); return([]));

m=n;
v=vector(l);
for(i=0, l-1, v[l-i]=m%d; m=floor(m/d););
return(v)
}

et demander à PARI/GP de lire ce fichier :
\r chiffres.gp

Il suffira ensuite d'appeler

Chiffres(3, 4, 28)

pour obtenir dans un vecteur de taille 4, la décomposition de 28 en base 3. Autrement dit,

[1, 0, 0, 1]

Finalement, ces exemples reflètent que peu la puissance de PARI/GP que je vous encourage à découvrir sur le site officiel, et qui est un formidable outil qui peut être très utile. Bien entendu vous pouvez calculer dans des corps finis, sur \mathbf Z/n\mathbf Z, manipuler des polynômes, des fractions rationnelles, les corps de nombres, les fonctions de Bessel, etc.)

Et pour finir, une jolie fonction tracée grâce à plot :

Le théorème de Wilson

Posted on by
7

Le théorème de Wilson est un résultat très connu d’artihmétique élémentaire, qui donne une caractérisation anecdotique des nombres premiers.

Dans cet article, je rappellerai les démonstrations connues et vous présenterai une démonstration un peu plus étonnante (et très anecdotique) basée sur le théorème de Sylow en théorie des groupes.

Historique

Le théorème de Wilson a été découvert par Ibn al-Haytham (aussi connu comme Alhazen) vers les années 1000, mais est attribué à John Wilson (un étudiant de Edward Waring) qui l’a conjecturé au dix-huitième sicècle. C’est Waring qui publia la conjecture 1770, mais sans démonstration. La première preuve est dûe à Lagrange, et a été publiée en 1773. Il semblerait que Leibniz connaissait le résultat un siècle plus tôt, mais ne l’a jamais publié.

Théorème de Wilson

Un entier naturel p>1 est premier si et seulement si

(p-1)!+1 \equiv 0 \pmod p

Cette caractérisation des nombres premiers est un peu anecdotique et ne constitue pas un test de primalité utile (calculer la factorielle de p-1 n’est pas envisageable). Son principal intérêt réside dans son histoire et la simplicité de son énoncé et de ses preuves (exceptée bien entendu celle dont je vais parler en dernier et qui justifie l’existence de cet article).

Démonstration du sens réciproque

Commençons par démontrer que si (p-1)!+1\equiv 0\pmod p, alors p est premier. Supposons que cette équation est vraie, cela signifie qu’il existe k\in\mathbf Z tel que

kp-(p-1)!=1

Ainsi, on a une relation de Bezout entre p et 1<\ell <p, donc p est premier.

Démonstration du sens direct

par Lagrange

Notons \mathbf F_p=\mathbf Z/p\mathbf Z. C’est un corps et le groupe \mathbf F_p^\times de ses inversibles est d’ordre p-1. Ainsi, le petit théorème de Fermat implique que ses p-1 éléments sont racines du polynôme X^{p-1}-1\in\mathbf F_p[X] de degré p-1. Comme \mathbf F_p est un corps, et que les éléments 1\leq \ell < p sont racines de ce polynôme, on a

X^{p-1}-1=(X-1)(X-2)\cdots(X-(p-1)).

En évaluant cette expression en p, on trouve le résultat voulu.

par Gauss

Si p=2, le résultat est clair. On suppose donc p impair.

Le principe de cette jolie démonstration est de grouper deux par deux les éléments de \mathbf F_p^\times. Pour tout k\in\mathbf F_p^\times, k\neq 1, -1, il existe \ell\in \mathbf F_p^\times, \ell\neq k (sinon k^2=1, i.e., k\in\left\{-1,1\right\} ce qu’on a justement exclu) tel que k\ell =1. Ainsi, lorsqu’on élimine, dans le produit \prod_{1\leq k\leq p-1}k\in\mathbf F_p, les paires d’inverses mutuels dont le produit vaut 1, il ne reste que

\prod_{1\leq k\leq p-1}k = 1\times (p-1) = -1 \in \mathbf F_p,

ce qui prouve le résultat.

… par le théorème de Sylow

Commençons par dénombrer le nombre de p-Sylow du groupe symétrique \mathfrak S_p

Dans \mathfrak S_p, il y a (p-1)! éléments d’ordre p. En effet, un élément d’ordre p est un p-cycle car p est premier, et un p-cycle peut s’écrire avec n’importe quel élément en premier, et c’est l’ordre des p-1 éléments restants qui importe.

Maintenant, un p-Sylow est de cardinal p^\alpha maximal avec p^\alpha qui divise le cardinal de \mathfrak S_p. Ainsi, \alpha=1 et les p-Sylow sont d’ordre p, et contiennent l’identité et p-1 p-cycles. Réciproquement, chaque p-cycle est contenu dans un p-Sylow. Il y a donc

\frac{(p-1)!}{p-1}=(p-2)! p-Sylow.

Avec le théorème de Sylow…

on sait que le nombre de p-Sylow est congru à 1 modulo p. Ainsi

(p-2)!\equiv 1\pmod p,

et en multipliant par p-1, on a

(p-1)!\equiv p-1\equiv -1\pmod p.

Petit mot de la fin

J’aime bien cette dernière démonstration pour deux raisons. La première, c’est qu’elle relie les sous-groupes de Sylow avec un théorème d’arithmétique élémentaire, ce qui prouve (si besoin est encore) qu’en mathématiques il n’y a jamais qu’un seul moyen de prouver les choses. La seconde est que les étudiants, face au dénombrement des p-Sylow de \mathfrak S_p, utilisent immédiatement le théorème de Sylow et restent bloqué avec

n_p\equiv 1\pmod p,

où n_p est le nombre de p-Sylow. Là, il faut remettre les mains dans le cambouis pour s’en sortir. De plus, ce n’est pas un exercice difficile.

L’historique et les premières démonstrations sont tirées de la page sur le théorème de Wilson de la wikipédia.

Simuler l’équation des ondes en dimension 2

Posted on by
1

On peut modéliser la propagation d’une onde dans un milieu homogène, linéaire et isotrope en dimension 2 par l’équation aux dérivées partielles suivante :

\dfrac{1}{c^2} \dfrac{\partial^2}{\partial t^2} u(t,x,y)=\dfrac{\partial^2}{\partial x^2} u(t,x,y)+\dfrac{\partial^2}{\partial y^2} u(t,x,y)

où u(t,x,y) représente l’amplitude de l’onde au temps t et au point de coordonnées (x,y) et c la vitesse de propagation de l’onde dans le milieu considéré. Nous allons voir comment simuler une telle équation en utilisant la méthode des différences finies. Dans la suite, on considère c=1, on se donne la condition initiale u^0(x,y)=u(0,x,y) et on suppose que l’onde n’a pas de vitesse initiale :

\dfrac{\partial }{\partial t} u(0,x,y)=0.

On fait ces hypothèses pour simplifier l’article, mais la méthode reste valable sans.

Regardons déjà ce qui se passe en dimension 1. Pour une fonction f régulière, la formule de Taylor nous assure que pour un réel \Delta x>0 petit on a

f(x+\Delta x)=f(x)+\Delta x f^\prime (x)+\dfrac{ \Delta x^2}{2} f^{\prime\prime }(x)+\mathcal O(\Delta x^3)

et

f(x-\Delta x)=f(x)-\Delta xf^\prime(x)+\dfrac{\Delta x^2}{2}f^{\prime\prime }(x)+\mathcal{O}(\Delta x^3).

En additionnant les deux lignes, on a donc

f(x+\Delta x)+f(x-\Delta x)=2f(x)+\Delta x^2f^{\prime\prime }(x)+\mathcal{O}(\Delta x^3).

On note f_i la valeur de f au point i\Delta x : f_i=f(i\Delta x). En négligeant le terme en \mathcal O(\Delta x^3) on a donc

f^{\prime\prime }_i =\dfrac{f_{i+1}+f_{i-1}-2f_i}{\Delta x^2}.

On peut faire le même raisonnement avec des fonctions de plusieurs variables. L’équation des ondes s’écrit alors

\dfrac{u_{i,j}^{n+1}+u_{i,j}^{n-1}-2u_{i,j}^n}{\Delta t^2}=\dfrac{u_{i+1,j}^n+u_{i-1,j}^{n}-2u_{i,j}^n}{\Delta x^2}+\dfrac{u_{i,j+1}^n+u_{i,j-1}^n-2u_{i,j}^n}{\Delta y^2}

où \Delta t, \Delta x et \Delta y sont des réels positifs petits et où u_{i,j}^n=u(n\Delta t,i\Delta x,j \Delta y). On peut écrire cette formule de manière différente :

u_{i,j}^{n+1}=2u_{i,j}^{n}-u_{i,j}^{n-1}+\dfrac{\Delta t^2}{\Delta x^2}(u_{i+1,j}^n+u_{i-1,j}^{n}-2u_{i,j}^n)+\dfrac{\Delta t^2}{\Delta y^2}(u_{i,j+1}^n+u_{i,j-1}^n-2u_{i,j}^n)

ce qui permet ainsi de calculer les valeurs approchées de u^{n+1} en fonction des valeurs de u^n et u^{n-1}. Il suffit donc d’écrire une boucle dans votre logiciel de calcul numérique favori (pour moi Matlab, sinon en gratuit il y a Scilab) pour calculer les solutions approchées. Il reste cependant deux problèmes. On se donne la condition de départ u^0. On ne peut calculer u^1 par la formule précédente, car il nous faut les expressions des deux termes u^0 et u^{-1} que l’on a pas! Heureusement Taylor nous vient encore une fois en aide, on a

u^1_{i,j}=u^0_{i,j}+\Delta t \dfrac{\partial}{\partial t} u^0_{i,j}+\dfrac{\Delta t^2}{2}\dfrac{\partial^2}{\partial t^2} u^0_{i,j}+\mathcal{O}(\Delta t^3).

La vitesse initiale est nulle, donc le terme d’ordre 1 est nul. De plus, on a

\dfrac{\partial^2}{\partial t^2} u^0_{i,j}=\dfrac{\partial^2}{\partial x^2} u^0_{i,j}+\dfrac{\partial^2}{\partial y^2} u^0_{i,j}

par l’équation des ondes. En refaisant encore une fois les approximations par Taylor, on trouve que

u_{i,j}^1=u_{i,j}^0+\dfrac{\Delta t^2}{2\Delta x^2}(u_{i+1,j}^0+u_{i-1,j}^{0}-2u_{i,j}^0)+\dfrac{\Delta t^2}{2\Delta y^2}(u_{i,j+1}^0+u_{i,j-1}^0-2u_{i,j}^0).

Le second problème est que l’ordinateur ne peut effectuer les calculs que sur un nombre fini de données. Il faut donc borner i et j. Cela pose un problème pour calculer la valeur de u aux bords du domaine. Si on prend -N \leq j \leq N, pour calculer u^n_{N,j} il faudrait connaitre u^{n-1}_{N+1,j}, ce que l’on ne connait pas. On se donne alors les valeurs sur le bord du domaine de calcul. On peut par exemple donner une valeur nulle aux bords, ce qui donnera une onde réfléchie sur les bords. On peut également développer des méthodes pour que les bords soient absorbants, bien qu’il subsiste toujours des réflexions parasites. Ici, j’ai évité le problème en arrêtant la simulation quand la première onde touche les bords.

Voilà ce que donne cette méthode implémentée sous Matlab. Ca ressemble pas si mal que ça à une mare dans laquelle on aurait jeté une pierre, non?

Vidéo de l’équation des ondes en dimension 2

Remarque : la vidéo chez moi est en rouge sous Windows Media Player et en bleu sous VLC, je n’ai aucune idée de la raison… Avez-vous une idée?

L’équation de transport unidimensionnelle

Posted on by
Reply

On appelle équation aux dérivées partielles une équation qui lie les dérivées partielles d’une ou plusieurs fonctions. Elles sont énormément utilisées en physique, dans des domaines aussi variées que la mécanique des fluides (équation de Navier-Stokes), l’électromagnétisme (équations de Maxwell), la thermodynamique (équation de la chaleur), la mécanique quantique (équation de Schrödinger), etc… On les trouve également dans des domaines comme la météorologie, la synthèse d’image ou l’économie  (équation de Black-Scholes).

On va étudier l’exemple le plus simple d’équation aux dérivées partielles : l’équation de transport unidimensionnelle. Elle s’écrit sous la forme

\displaystyle \frac{\partial}{\partial x} u(x,y)+c \frac{\partial}{\partial y} u(x,y)=0

où u(x,y) est la fonction inconnue des deux variables x et y et où c \neq 0 est une constante. Comme pour les équations différentielles ordinaires, il faut ajouter une condition initiale :

u(0,y)=u_0(y)

où u_0 est une fonction donnée de la variable y que l’on suppose assez régulière, par exemple dérivable à dérivée continue.

Résolution

On suppose que u est une solution régulière de notre problème, c’est-à-dire que u admet des dérivées partielles selon x et selon y et que ses dérivées partielles sont continues.

On va montrer que u est constante le long de certaines droites, appelées caractéristiques. Pour y_0 un réel, posons D(x)=cx+y_0 la droite de pente c qui passe par y_0.

On a par les règles de dérivation des fonctions composées :

\displaystyle\frac{\mathrm d}{\mathrm d x} u(x,D(x))=\frac{\partial}{\partial x} u(x,D(x))+D'(x) \frac{\partial}{\partial y}u(x,D(x)).

Or D'(x)=c, on a donc

\displaystyle\frac{\mathrm d}{\mathrm d x} u(x,D(x))=\frac{\partial}{\partial x} u(x,D(x))+c \frac{\partial}{\partial y}u(x,D(x))

et comme u est solution de l’équation de transport, on a

\displaystyle\frac{\mathrm d}{\mathrm d x} u(x,D(x))=0

ce qui montre que u est constante le long de toutes les caractéristiques. On a donc pour tout réel x

u(x,D(x))=u(0,D(0))=u(0,y_0)

ce qui s’écrit

u(x,cx+y_0)=u_0(y_0).

Soient x et y deux réels. Pour connaître la valeur de u(x,y), il suffit de connaitre la valeur de u_0 au point où la caractéristique passant par (x,y) coupe l’axe des ordonnées car u est constante le long de cette caractéristique. On cherche donc le réel y_0 tel que la caractéristique passant par (0,y_0) de pente c passe aussi par le point (x,y). Il faut donc que D(x)=cx+y_0, c’est-à-dire y_0=y-cx. Finalement, on obtient

\boxed{u(x,y)=u_0(y-cx).}

C’est la formule de Duhamel. Réciproquement, il est facile de vérifier que si on définit u par la formule de Duhamel, alors u admet des dérivées partielles continues et que

\begin{cases} \displaystyle\frac{\partial}{\partial x} u(x,y)+c \frac{\partial}{\partial y} u(x,y)=0 \\u(0,y)=u_0(y).\end{cases}.

On a donc montré que u définie par la formule de Duhamel est l’unique solution de l’équation de transport.

Un petit séjour à Singapour

Posted on by
1

Dans le cadre de mes études, j’ai effectué un stage à la Nanyang Technological University de Singapour en juin et juillet dernier, dans la School of Physical and Mathematical Science, c’est-à-dire le bâtiment réservé aux mathématiques, à la physique, à la chimie et à la chimie biologique. Cette université est très récente, et le bâtiment dans lequel je me trouvais était neuf (construit en 2007, 38000 m2, ouvert officiellement en juin 2009), très moderne, très grand, très beau : des conditions optimales pour travailler. Je n’étonnerai personne lorsque je révélerai que je n’ai rencontré qu’un seul Singapourien à l’université parmi les post-doc, professeurs, chercheurs et thésards. Cette université se sert du fait que Singapour est une pierre angulaire du commerce mondial pour en faire un rayonnement culturel international au niveau des sciences. Il est connu que les universités singapouriennes viennent aux congrès de mathématiciens (par exemple) avec des liste de 40 postes à pourvoir (ce qui est du jamais vu !).

SPMS, Singapore

J’ai rencontré de nombreuses personnes : des allemands, chinois, malais, indiens, turques, israéliens, suisses, américains, etc. C’est une expérience extraordinaire d’être en contact avec des gens de toute la terre, de discuter avec eux, de découvrir d’autres cultures, et de travailler ensemble. Des millions de dollars sont investis dans les universités à Singapour pour en faire des pôles d’excellence : la ville-État mise énormément sur la recherche et l’enseignement pour garder son avance sur les autres pays asiatiques, et s’immiscer au même niveau que l’Occident. Occident, qui soit-dit en passant, est très présent là bas dans la vie de tous les jours. Singapour est un étonnant mélange entre l’extrêmement moderne (Ion Orchard, Marina Bay Sands : le dernier hôtel casino, lubie du multi-milliardaire juif Sheldon Adelson) et le plus traditionnel, voire ressemblant aux pays en voie de développement (temples indiens, quartiers de petites maisons, superstition chinoise, etc.)… Mais le gouvernement lutte actuellement pour que les jeunes ne s’occidentalisent pas trop.

J’étais à l’université pendant les vacances scolaires de là bas, donc je n’ai pas vu l’effervescence des étudiants (excepté les remises de diplômes — correspondant à nos licences européennes — qui se déroulaient comme l’on voit dans les films américains), mais il y avait quand même beaucoup de personne dans l’immense campus, qui travaillent dans les différentes universités.

Au niveau des mathématiques, là-bas on travaille les mathématiques technologiques, pratiques. Dans l’autre université (très connue) de Singapour, la NUS, il y a plus de mathématiques « pures ». Le département de cryptologie est très développé et on trouve de nombreuses personnes qui travaillent dessus : mathématiciens, informaticiens, ingénieurs, etc. Cependant malgré cette effervescence, d’un point de vue extérieur, le cadre ne semblait pas du tout stressant : les conditions de travail sont très bonnes — bureaux clairs, environnement vert, climatisation (ô combien nécessaire), une grande bibliothèque, proximité des laboratoires de recherche avec les bureaux personnels des professeurs, salles de discussions, nombreux séminaires, et un calme qui vous transperce complètement et propice au travail…

Là-bas j’ai travaillé sur le partage de secret en cryptologie, qui a fait l’objet d’un article précédent sur ce blog.

Un séjour très instructif et agréable (excepté le climat pas très accueillant : chaud et humide, très humide, trop humide… équatorial !).