Dans cet article, on va proposer une réponse à la question posée dans l’article intitulé « Un arbre pythagoricien » par Etienne Ghys et Jos Leys, publié dans la revue Images des mathématiques, CNRS, 2013.
Cet article peut être librement diffusé à l’identique dans la limite d’une utilisation non commerciale suivant la licence CC-nc-nd
(http://creativecommons.org/licenses/by-nc-nd/3.0/fr/)
Le problème à résoudre
Résumons le propos de l’article (dont on a repris sans vergogne les figures) :
L’arbre de Pythagore est une figure fractale que l’on construit assez naturellement par récursivité.
Pour l’arbre de Pythagore classique en 2D, le scénario de construction est le suivant :
On part d’un carré et on en constuit deux autres...
puis quatre autres...
et on recommence...
Les branches de l’arbre de Pythagore viennent vite à se chevaucher.
L’article démontre mathématiquement ce chevauchement, que l’on constate de visu.
On passe alors à un premier arbre de Pythagore 3D en remplaçant les carrés par des cubes que l’on dispose naturellement en faisant en sorte que le plan frontal soit un plan de symétrie.
Enfin, on opére une rotation de 90° à chaque étape pour obtenir un arbre plus intéressant :
Visuellement, on constate que pour cet arbre de Pythagore 3D il n’y a plus chevauchement des branches.
Le problème posé est de le démontrer.
La représentation par CaRMetal et DGPad
Pour représenter cet arbre de Pythagore 3D, Etienne Ghys et Jos Leys ont utilisé POV-RAY, qui génère de très belles images 3D statiques.
On va ici représenter l’arbre véritablement en 3D (dans un repère mobile) en utilisant les logiciels CaRMetal et DGPad.
Ces deux logiciels sont développés par Eric Hakenholz, et leur gestion de la 3D est différente.
Dans le cas présent, CaRMetal est plus robuste et permet d’explorer des arbres plus profonds.
DGPad quant à lui est remarquable par son approche radicalement différente de la 3D et a le très grand intérêt de pouvoir aisément insérer en ligne des figures dynamiques.
Les copies d’écran seront moins convaincantes esthétiquement (il n’y aura pas d’effet d’ombrage), mais en construisant l’arbre de cette façon on pourra accéder aux composantes numériques de l’arbre.
Cette représentation est plus que visuelle, elle donnera le squelette de l’arbre dont on n’avait que l’ADN.
Commençons par décrire la gestion de la 3D dans CaRMetal.
La 3D dans CaRMetal
CaRMetal est un logiciel de pseudo-3D : en mode 3D, on dispose d’un trièdre mobile (O,X,Y,Z) pilotable par clic droit dans l’espace 2D.
En projetant sur ce trièdre, on peut représenter n’importe quel point 3D.
Cette approche peut alors tirer parti de la puissance et de la simplicité des scripts CaRMetal.
Ces scripts comprennent :
- les CaRScripts qui traduisent les outils de CaRMetal ;
- une bonne partie du langage JavaScript.
Le CaRScript Point(s,t) crée un point de coordonnées (s,t) par effet de bord.
s et t peuvent contenir des expressions CaRMetal. Par exemple, Point(« x(A) », « y(A) ») construit un doublon du point A (qui doit déjà être construit pour que les expressions x(A) et y(A) soient reconnues).
Pour un point M(a,b,c) de l’espace, on a : $M = O + a\overrightarrow{OX} +b\overrightarrow{OY} + c\overrightarrow{OZ}$.
On peut alors créer le projeté d’un point de l’espace de coordonnées (a,b,c) par l’instruction suivante :
Point("x(O)+("+a+")*(x(X)-x(O))+("+b+")*(x(Y)-x(O))+("+c+")*(x(Z)-x(O))","y(O)+("+a+")*(y(X)-y(O))+("+b+")*(y(Y)-y(O))+("+c+")*(y(Z)-y(O))");
Script de construction de l’arbre 3D
Pour simplifier l’écriture du script, on va mettre en place des outils réutilisables.
Le premier outil est une fonction (dans l’esprit, un CaRScript) Point3D.
Point3D
- var Point3D = function (a,b,c) {
- var that= {};
- that.x = a;
- that.y = b;
- that.z =c;
- that.name= Point("x(O)+("+a+")*(x(X)-x(O))+("+b+")*(x(Y)-x(O))+("+c+")*(x(Z)-x(O))","y(O)+("+a+")*(y(X)-y(O))+("+b+")*(y(Y)-y(O))+("+c+")*(y(Z)-y(O))");
- that.getX3D = function () {
- return that.x;
- }
- that.getY3D = function () {
- return that.y;
- }
- that.getZ3D = function () {
- return that.z;
- }
- that.getName2D = function () {
- return that.name;
- }
- return that;
- };
Cette fonction s’utilise comme le CaRScript Point() mais crée un point 3D.
La propriété name donne le nom du projeté dans l’espace 2D.
Dans certains cas (=pour les calculs), on ne voudra pas que le point 2D soit créé par effet de bord. On va donc créer une fonction Point3Dc analogue (convention : suffixe c pour calcul), mais sans création du projeté :
- var Point3Dc = function (a,b,c) {
- var that= {};
- that.x = a;
- that.y = b;
- that.z =c;
- that.getX3D = function () {
- return that.x;
- }
- that.getY3D = function () {
- return that.y;
- }
- that.getZ3D = function () {
- return that.z;
- }
- return that;
- };
On va ensuite créer une fonction somme sur les points 3D, avec ou sans création du projeté.
Somme
- var somme = function(a,b) {
- return Point3D(a.getX3D()+b.getX3D(),a.getY3D()+b.getY3D(),a.getZ3D()+b.getZ3D());
- };
- var sommec = function(a,b) {
- return Point3Dc(a.getX3D()+b.getX3D(),a.getY3D()+b.getY3D(),a.getZ3D()+b.getZ3D());
- };
On crée de même une fonction vect, un produit par un réel, et un produit vectoriel :
Vecteur
- var vect = function(a,b) {
- return Point3D(b.getX3D()-a.getX3D(),b.getY3D()-a.getY3D(),b.getZ3D()-a.getZ3D());
- };
- var vectc = function(a,b) {
- return Point3Dc(b.getX3D()-a.getX3D(),b.getY3D()-a.getY3D(),b.getZ3D()-a.getZ3D());
- };
Produit par un réel
- var prodR = function(a,b) {
- return Point3D(a*b.getX3D(),a*b.getY3D(),a*b.getZ3D());
- };
- var prodRc = function(a,b) {
- return Point3Dc(a*b.getX3D(),a*b.getY3D(),a*b.getZ3D());
- };
Produit vectoriel
- var prodVect = function(e,g) {
- var nx=e.getY3D()*g.getZ3D()- e.getZ3D()*g.getY3D();
- var ny=e.getZ3D()*g.getX3D()- e.getX3D()*g.getZ3D();
- var nz=e.getX3D()*g.getY3D()- e.getY3D()*g.getX3D();
- return Point3D(nx,ny,nz);
- };
- var prodVectc = function(e,g) {
- var nx=e.getY3D()*g.getZ3D()- e.getZ3D()*g.getY3D();
- var ny=e.getZ3D()*g.getX3D()- e.getX3D()*g.getZ3D();
- var nz=e.getX3D()*g.getY3D()- e.getY3D()*g.getX3D();
- return Point3Dc(nx,ny,nz);
- };
Enfin, on construit une fonction norme :
Norme
- var norme = function(u,v) {
- return Math.sqrt(vectc(u,v).getX3D()*vectc(u,v).getX3D()+vectc(u,v).getY3D()*vectc(u,v).getY3D()+vectc(u,v).getZ3D()*vectc(u,v).getZ3D());
- };
On crée ensuite une fonction cube qui s’applique à trois points e, f, g de l’espace 3D
(les points doivent former un triangle rectangle isocèle pour donner un cube).
On utilise les outils créés précédemment.
- var cube = function (e,f,g) {
- var that= {};
- that.s1 = e;
- that.s2= f;
- that.s3 =g;
- that.s4= vect(f,sommec(e,g));
- var areteNormale= prodVectc(vectc(f,e),vectc(f,g));
- var a= norme(e,f);
- areteNormale= prodRc(1/a,areteNormale);
- that.s5= somme(e,areteNormale);
- that.s6= somme(f,areteNormale);
- that.s7= somme(g,areteNormale);
- that.s8= somme(that.s4,areteNormale);
- that.getS4 = function () {
- return that.s4;
- }
- that.getS5 = function () {
- return that.s5;
- }
- that.getS6 = function () {
- return that.s6;
- }
- that.getS7 = function () {
- return that.s7;
- }
- that.getS8 = function () {
- return that.s8;
- }
- Polygon(that.s1.name+","+that.s2.name+","+that.s3.name+","+that.s4.name);
- Polygon(that.s1.name+","+that.s2.name+","+that.s6.name+","+that.s5.name);
- Polygon(that.s2.name+","+that.s3.name+","+that.s7.name+","+that.s6.name);
- Polygon(that.s3.name+","+that.s4.name+","+that.s8.name+","+that.s7.name);
- Polygon(that.s4.name+","+that.s1.name+","+that.s5.name+","+that.s8.name);
- Polygon(that.s5.name+","+that.s6.name+","+that.s7.name+","+that.s8.name);
- return that;
- };
On s’intéresse ensuite à la récursion.
Voici le schéma (n pour nouveau) :
La fonction récursive est la suivante :
- function arbre(o,p,q,n) {
- if (n>0) {
- var c= cube(o,p,q);
- var no1=somme(prodRc(0.5+0.5/Math.sqrt(2),c.getS8()),prodRc(0.5-0.5/Math.sqrt(2),c.getS7()));
- var inter1=somme(prodRc(0.5+0.5/Math.sqrt(2),c.getS7()),prodRc(0.5-0.5/Math.sqrt(2),c.getS8()));
- var np1= somme(no1,sommec(prodRc(0.5,vectc(c.getS7(),c.getS6())),prodRc(0.5,vectc(o,c.getS5()))));
- var nq1 = somme(inter1,sommec(prodRc(0.5,vectc(c.getS7(),c.getS6())),prodRc(0.5,vectc(o,c.getS5()))));
- var no2 = somme(prodRc(0.5+0.5/Math.sqrt(2),c.getS6()),prodRc(0.5-0.5/Math.sqrt(2),c.getS5()));
- arbre(no1,np1,nq1,n-1);
- arbre(no2,nq1,np1,n-1);
- }
- }
Remarques :
- inter1 est un point intermédiaire sur l’arête, qui permet de construire nq1 ;
- inter1 et no1 sont obtenus par barycentre des points S7 et S8 (coef $\dfrac{1}{2} \pm \dfrac{1}{2\sqrt{2}}$).
On termine en lançant la fonction arbre :
- niv=8;
- arbre(Point3D(-1,1,-1),Point3D(1,1,-1),Point3D(1,-1,-1),niv);
Et on obtient bien l’arbre (choisir une forme « point »... pour les points) :
Avant d’aller plus loin, on va procéder à quelques modifications que m’ont inspirées mes fidèles relecteurs que sont Yves Martin et Alain Busser :
* On ne va pas tracer les points par défaut. Il suffit pour cela d’ajouter cette ligne de code dans la fonction Point3D :
- SetHide(that.name,true);
On obtient un arbre qui n’est plus « pollué » par les « moucherons »...
* On va s’intéresser à la lisibilité des scripts qui traduisent des égalités mathématiques de géométrie 3D.
Comment faire sachant que l’on ne peut pas surcharger les opérateurs en Javascript ?...
On va opter pour une notation polonaise inversée.
La fonction calcul3D utilisera les fonctions précédentes en notation polonaise inversée.
Voici le code de cette fonction :
- var calcul3DTab = function(a) {
- var b=[];
- if (a.length==3){
- return a[2](a[0],a[1]);
- }
- else {
- var i=0;
- do {
- b[i]=a[i];
- i++;
- }
- while ((typeof a[i])!="function");
- b[i-2]= a[i](a[i-2],a[i-1]);
- i--;
- while (i<a.length-2) {
- b[i]= a[i+2];
- i++;
- }
- return calcul3DTab(b);
- }
- };
- var calcul3D = function() {
- return calcul3DTab(arguments);
- };
Le code des fonctions cube et arbre devient alors :
- var cube = function (e,f,g) { // e,f,g sommets coplanaires "consécutifs"
- var that= {};
- that.s1 = e;
- that.s2= f;
- that.s3 =g;
- that.s4= calcul3D(f,e,g,sommec,vect);
- var areteNormale= calcul3D(f,e,vectc,f,g,vectc,prodVectc);
- var a= norme(e,f);
- areteNormale= prodRc(1/a,areteNormale);
- that.s5= somme(e,areteNormale);
- that.s6= somme(f,areteNormale);
- that.s7= somme(g,areteNormale);
- that.s8= somme(that.s4,areteNormale);
- that.getS4 = function () {
- return that.s4;
- }
- that.getS5 = function () {
- return that.s5;
- }
- that.getS6 = function () {
- return that.s6;
- }
- that.getS7 = function () {
- return that.s7;
- }
- that.getS8 = function () {
- return that.s8;
- }
- Polygon(that.s1.name+","+that.s2.name+","+that.s3.name+","+that.s4.name);
- Polygon(that.s1.name+","+that.s2.name+","+that.s6.name+","+that.s5.name);
- Polygon(that.s2.name+","+that.s3.name+","+that.s7.name+","+that.s6.name);
- Polygon(that.s3.name+","+that.s4.name+","+that.s8.name+","+that.s7.name);
- Polygon(that.s4.name+","+that.s1.name+","+that.s5.name+","+that.s8.name);
- Polygon(that.s5.name+","+that.s6.name+","+that.s7.name+","+that.s8.name);
- return that;
- };
- function arbre(o,p,q,n) {
- if (n>0) {
- var c= cube(o,p,q);
- var no1= calcul3D(0.5+0.5/Math.sqrt(2),c.getS8(),prodRc,0.5-0.5/Math.sqrt(2),c.getS7(),prodRc,somme);
- var inter1=calcul3D(0.5+0.5/Math.sqrt(2),c.getS7(),prodRc,0.5-0.5/Math.sqrt(2),c.getS8(),prodRc,somme);
- var np1= calcul3D(no1,0.5,c.getS7(),c.getS6(),vectc,prodRc,0.5,o,c.getS5(),vectc,prodRc,sommec,somme);
- var nq1= calcul3D(inter1,0.5,c.getS7(),c.getS6(),vectc,prodRc,0.5,o,c.getS5(),vectc,prodRc,sommec,somme);
- var no2= calcul3D(0.5+0.5/Math.sqrt(2),c.getS6(),prodRc,0.5-0.5/Math.sqrt(2),c.getS5(),prodRc,somme);
- arbre(no1,np1,nq1,n-1);
- arbre(no2,nq1,np1,n-1);
- }
- }
La 3D dans DGPad
Avant d’aller plus loin, on va donner une figure dynamique manipulable en ligne réalisée avec DGPad (arbre de niveau 5).
(Clic droit sur ordinateur ou 2 doigts sur tablette pour faire tourner le repère)
Voici pour information le script de cette figure :
Le script de CaRMetal a été adapté pour DGPad.
On utilise des expressions pour générer les points 3D. Et la transposition des calculs géométriques au niveau des expressions permet en quelque sorte de surcharger les opérateurs. Il faut ensuite manipuler les expressions pour parvenir au résultat attendu.
On n’en dira pas davantage ici sur la 3D dans DGPad. Vous trouverez plus d’information sur les expressions de DGPad dans cet article du même numéro et plus spécifiquement pour la 3D dans cette vidéo Youtube d’Eric Hakenholz et dans cet article détaillé d’Yves Martin sur le site de l’IREM de la Réunion.
Démonstration qu’il n’y a pas chevauchement
Maintenant que l’on dispose de la fonction arbre, on va pouvoir démontrer qu’il n’y a pas chevauchement.
On va partir des mêmes points, donc avec une arête L=2.
La question est de savoir si l’on traverse le plan y=0.
Pour les sommets déjà tracés, on peut aisément demander à CaRMetal si c’est le cas en lui faisant calculer (à l’intérieur de la fonction arbre) la valeur minimale de abs($y_M$) pour tout sommet M d’un cube.
Pour les sommets qui ne sont pas tracés, on utilise le fait (démontré dans l’article cité) que ces sommets sont inclus dans une boule de centre le centre du cube bourgeon et de rayon $ R \approx 3,82L \approx 4L$.
Si cette boule ne traverse pas le plan y=0, il en sera de même pour tous les futurs sommets.
Pour visualiser la situation, on va même modifier le script pour faire afficher ces boules en perspective, leur projeté n’étant pas trop difficile à tracer dans l’espace 2D (ce sont des disques de même rayon).
Voici la nouvelle fonction arbre :
- function arbre(o,p,q,n) {
- if (n>0) {
- var c= cube(o,p,q);
- var no1= calcul3D(0.5+0.5/Math.sqrt(2),c.getS8(),prodRc,0.5-0.5/Math.sqrt(2),c.getS7(),prodRc,somme);
- var inter1=calcul3D(0.5+0.5/Math.sqrt(2),c.getS7(),prodRc,0.5-0.5/Math.sqrt(2),c.getS8(),prodRc,somme);
- var np1= calcul3D(no1,0.5,c.getS7(),c.getS6(),vectc,prodRc,0.5,o,c.getS5(),vectc,prodRc,sommec,somme);
- var nq1= calcul3D(inter1,0.5,c.getS7(),c.getS6(),vectc,prodRc,0.5,o,c.getS5(),vectc,prodRc,sommec,somme);
- var no2= calcul3D(0.5+0.5/Math.sqrt(2),c.getS6(),prodRc,0.5-0.5/Math.sqrt(2),c.getS5(),prodRc,somme);
- if (n==1) {
- var cDisque= prodR(0.5,sommec(o,c.getS7()));
- var disque= FixedCircle(cDisque.name,8/Math.pow(Math.sqrt(2),niv-1));
- SetFilled(disque,true);
- }
- arbre(no1,np1,nq1,n-1);
- arbre(no2,nq1,np1,n-1);
- }
- }
Et le résultat pour 5 niveaux :
Voici le résultat avec DGPad :
On voit que la méthode appiquée au niveau 5 n’est pas probante, les boules se chevauchant.
Mais dès le niveau 6, on voit que la méthode va fonctionner (c’est « limite », mais en zoomant, on voit que les boules ne se chevauchent pas).
Il ne reste plus qu’à le vérifier par le calcul.
On crée une variable globale distanceMini.
Cette variable contiendra la valeur minimale de abs($y_M$) [1] pour tous les sommets M de cubes (on enlève les deux premiers) et de $\text{abs}(y_C) –$ rayonBoule pour les centres C d’un cube terminal.
Si distanceMini reste supérieur à 0 pour un arbre de niveau donné (quel qu’il soit), on est assuré qu’il n’y a pas chevauchement pour l’arbre complet.
Voici la dernière fonction arbre :
- function arbre(o,p,q,n) {
- if (n>0) {
- var c= cube(o,p,q);
- var no1= calcul3D(0.5+0.5/Math.sqrt(2),c.getS8(),prodRc,0.5-0.5/Math.sqrt(2),c.getS7(),prodRc,somme);
- var inter1=calcul3D(0.5+0.5/Math.sqrt(2),c.getS7(),prodRc,0.5-0.5/Math.sqrt(2),c.getS8(),prodRc,somme);
- var np1= calcul3D(no1,0.5,c.getS7(),c.getS6(),vectc,prodRc,0.5,o,c.getS5(),vectc,prodRc,sommec,somme);
- var nq1= calcul3D(inter1,0.5,c.getS7(),c.getS6(),vectc,prodRc,0.5,o,c.getS5(),vectc,prodRc,sommec,somme);
- var no2= calcul3D(0.5+0.5/Math.sqrt(2),c.getS6(),prodRc,0.5-0.5/Math.sqrt(2),c.getS5(),prodRc,somme);
- if (n==1) {
- var cDisque= prodR(0.5,sommec(o,c.getS7()));
- var disque= FixedCircle(cDisque.name,8/Math.pow(Math.sqrt(2),niv-1));
- SetFilled(disque,true);
- if (Math.abs(cDisque.getY3D())-8/Math.pow(Math.sqrt(2),niv-1)<distanceMini) {
- distanceMini=Math.abs(cDisque.getY3D())-8/Math.pow(Math.sqrt(2),niv-1);
- }
- }
- else {
- if (n<niv-1) {
- if (Math.abs(o.getY3D())<distanceMini) distanceMini=Math.abs(o.getY3D());
- if (Math.abs(p.getY3D())<distanceMini) distanceMini=Math.abs(p.getY3D());
- if (Math.abs(q.getY3D())<distanceMini) distanceMini=Math.abs(q.getY3D());
- if (Math.abs(c.getS4().getY3D())<distanceMini) distanceMini=Math.abs(c.getS4().getY3D());
- if (Math.abs(c.getS5().getY3D())<distanceMini) distanceMini=Math.abs(c.getS5().getY3D());
- if (Math.abs(c.getS6().getY3D())<distanceMini) distanceMini=Math.abs(c.getS6().getY3D());
- if (Math.abs(c.getS7().getY3D())<distanceMini) distanceMini=Math.abs(c.getS7().getY3D());
- if (Math.abs(c.getS8().getY3D())<distanceMini) distanceMini=Math.abs(c.getS8().getY3D());
- }
- }
- arbre(no1,np1,nq1,n-1);
- arbre(no2,nq1,np1,n-1);
- }
- else Println("distance minimale au plan y=0 : "+ distanceMini);
- }
Cela étant, il y a une petite arnaque dans le code précédent, c’est l’utilisation de la fonction Math.abs (il y a peu de chance que l’on obtienne une valeur négative en utilisant cette fonction…)
En fait, il faut remplacer cette fonction par Id quand on est sur la branche de droite et -Id quand on est sur la branche de gauche.
On modifie le code en conséquence en lançant d’abord une fonction départ qui permet de distinguer les deux branches principales. Et on ajoute un paramètre abs$=\pm1$ à la fonction arbre qui permet de savoir sur quelle branche on se trouve.
Le résultat du script (donné dans la dernière figure du fichier CaRMetal joint) est le même.
Pour 6 niveaux, on obtient une « distance minimale » de 0,00996.
Par conséquent les arbres ne se chevauchent pas.
On pourra légitimement préférer une démonstration mathématique, mais d’un certain point de vue il s’agit bien d’une preuve.
On pourra noter que l’argument mathématique de la boule, qui permettait de démontrer le chevauchement pour l’arbre 2D, a été greffé sur l’arbre 3D pour montrer le non chevauchement.
On joint un fichier CaRMetal avec quatre figures (= avec quatre scripts) :
- arbre
- arbre + boules
- arbre+boules+distance mini
- arbre+boules+distance mini sans arnaque
On choisit la figure par onglet, puis il suffit de lancer le script arbre.