retour index

chapitre 9

Oscillateurs chaotiques

Pendule forcé amorti

Oscillateur anharmonique forcé

Le pendule paramétrique

Aiguille dans deux champs magnétiques

Trajectoires d'une boule de billard

 

Une autre famille de systèmes pouvant présenter des régimes chaotiques est constituée par les oscillateurs. On se souvient que, pour avoir un mouvement chaotique, il faut une application non linéaire et un espace de phase de dimension au moins égale à 3. Pour les oscillateurs étudiés habituellement à un seul degré de liberté de position, il faudra un mouvement forçé conduisant à une équation différentielle du second ordre avec un second membre dépendant du temps, équation non autonome que l'on transforme en un système autonome à 3 dimensions, donc susceptible de donner du chaos.

Pendule amorti forcé

On sait que l'équation du pendule simple constitué d'une masse au bout d'un fil accroché en un point fixe et sans frottement est de la forme
MATH

l est la longueur du pendule et g l'accélération de la pesanteur. C'est un mouvement à 2 degrés de liberté (il faut 2 conditions initiales pour résoudre l'équation). Le mouvement d'un pendule simple est toujours régulier. S'il y a des frottements, le mouvement s'amortit et la masse retourne à sa position d'équilibre, qui est un point fixe attracteur. On pourra par contre avoir un système dissipatif pouvant donner des régimes chaotiques si l'on ajoute un frottement et un entretien ;l'équation sera alors
MATH

où 2$\lambda $ est le coefficient d'amortissement et MATH est la pulsation propre du système. On posera $\theta =2\pi x$ pour qu'un tour complet corresponde à x = 1, et $\omega =2\pi $ de sorte que l'entretien ait une période unité ;on a alors MATH et MATH

L'équation s'écrit :
MATH

Tracé de l'attracteur étrange

x étant (à un facteur près) l'angle du pendule avec la verticale, on peut identifier x avec x + 1 qui représente le même point. On prendra en fait la partie fractionnaire de x si x > 0 et la partie fractionnaire de x augmentée de 1 si x < 0 de façon à identifier les frontières gauche et droite du tracé comme si celui-ci était tracé sur un cylindre de périmètre 1.

On trace l'attracteur de la carte de période-1, qui est une stroboscopie de période égale à celle de l'entretien, c'est-à-dire que l'on va considérer les points de l'espace de phase atteints lorsque le temps est un multiple de la période, par exemple pour t = 0, 1, 2, 3.

On prend l'équation sous la forme MATH qui est bien de la forme précédente (on pose $2\lambda =2\pi c,$ $\omega _{0}=2\pi $ et MATH).

On choisit c = 0,2 et $\rho =2,5.$

> restart: with(DEtools):with(plots):

> c:=0.2:rho:=2.5:p:=phaseportrait([D(x)(t)=u,

D(u)(t)=6.28*(-c*u-sin(2*Pi*x)+rho*sin(2*Pi*t))],[x,u],

t=0..3000,[[x(0)=1,u(0)=0]],stepsize=0.1,

linecolor=black,scene=[x,u]):

> a:=op(1,op(1,p)):(on extrait la liste des points)

> liste1:=op([]):for i to (nops(a)-1)/10 do

if op(1,op(10*i,a))>0

then liste1:=liste1,[frac(op(1,op(10*i,a))),op(2,op(10*i,a))]

else liste1:=liste1,[1+frac(op(1,op(10*i,a))),

op(2,op(10*i,a))] fi od:(on prend la partie fractionnaire si x>0, sinon 1+la partie fractionnaire)

> plot([liste1],style=point):


chap9__17.pngfigure 9.1 : attracteur étrange de la carte de temps-1 du pendule amorti forcé avec $c=0.2$ et $\rho =2,5$

La structure fractale en couches fines de l'attracteur se voit d'autant mieux que le nombre de points est plus grand, ce qui allonge d'autant plus le temps de calcul. Cette structure est caractéristique du phénomène d'étirement-repliement, opération de mélange menant au chaos.

Comportement en fonction de l'amplitude du forçage

Revenons à une amplitude de forçage faible ;après un régime transitoire, le pendule oscille régulièrement de part et d'autre de sa position d'équilibre, ce qui se traduit par une ellipse dans le plan de phase (x, u). Si le forçage augmente, l'amplitude des oscillations augmente aussi, et, si celle-ci dépasse un demi-tour, le pendule va pouvoir effectuer des tours complets :le pendule décroche.

On observe alors un régime transitoire chaotique, dont le pendule sort pour se fixer dans un des trois régimes possibles qui sont soit un régime d'oscillations régulières, soit un régime où le pendule tourne régulièrement dans le même sens à la fréquence du forçage, dans un sens ou dans l'autre. Le choix du régime final dépend des conditions initiales ;il y a trois bassins d'attraction :

> restart: with(DEtools):

> c:=0.2:rho:=1.70:phaseportrait([D(x)(t)=u,

D(u)(t)=2*Pi*(-c*u-sin(2*Pi*x)+rho*sin(2*Pi*t))],[x,u],t=0..8,

[[x(0)=0.1,u(0)=0],[x(0)=-0.1,u(0)=0],[x(0)=0.5,u(0)=-1]],

stepsize=0.1,linecolor=red,scene=[t,x]);(3 conditions initiales)


chap9__20.pngfigure 9.2 : $c=0,2$ ; $\rho =1,7$ ; on a 3 possibilités d'état final différent

On a 3 régimes possibles pour les mêmes valeurs des paramètres selon l'état initial et donc trois attracteurs. Il y a multistabilité.

Le comportement est en fait complexe et très sensible aux valeurs des paramètres. Pour $\rho $ voisin de 2,17, le système n'arrive pas à se fixer sur un des trois régimes précédents et l'attracteur devient chaotique suivant un scénario qui rappelle une intermittence :

> restart: with(DEtools):

> c:=0.2:rho:=2.171:phaseportrait([D(x)(t)=u,

D(u)(t)=2*Pi*(-c*u-sin(2*Pi*x)+rho*sin(2*Pi*t))],[x,u],

t=100..800,[[x(0)=0.5,u(0)=0]],stepsize=0.1,

linecolor=red,scene=[t,u]);


chap9__23.pngfigure 9.3 : pour $\rho =2,171 $ , on a un régime intermittent

Avec l'option scene = [t, x], on obtient :


chap9__24.pngfigure 9.4 : régime intermittent :x en fonction de t

Le système tourne un certain temps dans un sens,puis change brutalement de sens.

Traçons l'attracteur de la carte de période-1 par le même programme que celui utilisé pour l'attracteur étrange, le temps étant compris entre 100 et 800 :


chap9__25.pngfigure 9.5 :attracteur du régime intermittent

Il faut comparer cette figure avec la figure 9.1 dont elle est l'esquisse. On voit apparaitre le dessin de l'attracteur étrange avec sa structure interne.

Si on augmente le forçage, les phases chaotiques se resserrent et le régime devient franchement chaotique.

Pour des valeurs supérieures du forçage, on retrouve des régimes réguliers, puis à nouveau chaotiques pour des valeurs encore supérieures.

début chapitre

Oscillateur anharmonique forcé

Oscillateur libre

C'est l'oscillateur à deux puits, appelé encore oscillateur de Duffing.

L'équation de l'oscillateur anharmonique libre non amorti est de la forme
MATH

En multipliant par $\overset{.}{x}$ et en intégrant, il vient
MATH

qui s'interprète comme la conservation de l'énergie mécanique $E, $ somme de l'énergie cinétique MATH(la masse étant prise égale à 1) et de l'énergie potentielle MATH, représentée comme suit :

>plot(x^4/4-x^2/2,x=-1.6..1.6);


chap9__32.pngfigure 9.6 : énergie potentielle de l'oscillateur anharmonique

Il s'agit bien d'un oscillateur à deux puits :on a deux points fixes pour $x=\pm 1$ qui sont des centres pour l'oscillateur non amorti et des puits si celui-ci est amorti, et un point-selle pour x = 0.

Traçons le portrait de phase en prenant plusieurs conditions initiales :

> restart: with(DEtools):

> phaseportrait([D(x)(t)=u,D(u)(t)=x-x^3],[x,u],t=0..12,

[seq([x(0)=1,u(0)=i/10],i=[1,3,5,7,9,11]),

seq([x(0)=-1,u(0)=i/10],i=[1,3,5,7,8])],stepsize=0.1,

scene=[x,u],linecolor=red,arrows=none);


chap9__34.pngfigure 9.7 : portrait de phase de l'oscillateur anharmonique libre non amorti

Considérons maintenant l'oscillateur amorti, d'équation


MATH

Traçons le portrait de phase :

>phaseportrait([D(x)(t)=u,D(u)(t)=x-x^3-0.1*u],[x,u],t=0..60,

[[x(0)=2,u(0)=0]],stepsize=0.1,scene=[x,u],

linecolor=red,arrows=none);


chap9__36.pngfigure 9.8 : portrait de phase de l'oscillateur anharmonique amorti

Du fait des frottements, l'énergie mécanique diminue, le point tourne dans le sens des aiguilles en se rapprochant ;à un instant donné, il entre dans l'un des deux puits et tourne autour de celui-ci jusqu'à sa position finale au fond du puits. Le système dépendant de deux paramètres (deux degrés de liberté), il n'y a pas de chaos possible.

Oscillateur forcé

Ce n'est plus le cas si l'on introduit un forçage sinuoïdal, comme c'est le cas pour l'équation


MATH

Traçons le nouveau portrait de phase,avec $\rho =0,3$ :

>phaseportrait([D(x)(t)=u,D(u)(t)=x-x^3-0.1*u+0.3*sin(t)],

[x,u],t=0..60,[[x(0)=1.8,u(0)=0]],stepsize=0.1,scene=[x,u],

linecolor=red,arrows=none);


chap9__40.pngfigure 9.9 : portrait de phase de l'équation :MATH

On voit que la trajectoire est chaotique. L'oscillateur tourne tantôt autour d'un point fixe, tantôt autour de l'autre ou autour des deux sans parvenir à se fixer.

On peut tracer l'attracteur de la carte de période-1 de l'oscillateur. Il faut pour cela prendre un forçage de période 1 de la forme $\rho \sin (2\pi t)$ ;cela revient à diviser MATH par 4$\pi ^{2}$ et $\overset{.}{x}=u$ par 2$\pi $ ;on prend $\rho =3$:

> restart: with(DEtools):with(plots):

> p:=phaseportrait([D(x)(t)=u,

D(u)(t)=4*Pi^2*(x-x^3-0.1*u/(2*Pi)+3*sin(2*Pi*t))],[x,u],

t=0..1000,[[x(0)=1.8,u(0)=0]],stepsize=0.05,scene=[x,u],

linecolor=red):

> a:=op(1,op(1,p)):nops(a);

> plot([seq(op(20*i,a),i=1..(nops(a)-1)/20)],style=point);


chap9__47.pngfigure 9.10 : attracteur étrange de Duffing

L'attracteur s'épaissit si l'on diminue le frottement.

On peut également tracer en une seule fois 4 cartes de temps-1 à des instants différents, en traçant avec pointplot3d un point sur 5, correspondant à 4 points par période puisque le pas d'intégration est de 0,05 et la période de 1. Les points de la même carte sont tous mis à la même abcisse en prenant la partie fractionnaire du temps :

> restart:with(DEtools):with(plots):

> p:=DEplot3d([D(x)(t)=u,

D(u)(t)=4*Pi^2*(x-x^3-0.1*u/(2*Pi)+3*sin(2*Pi*t))],[x,u],

t=0..500,[[x(0)=1.8,u(0)=0]],stepsize=0.05,scene=[t,x,u],

linecolor=red):

> a:=op(1,op(1,p)):nops(a);

> pointplot3d([seq([frac(op(1,op(5*i,a))),op(2,op(5*i,a)),

op(3,op(5*i,a))],i=1..(nops(a)-1)/5)]);


chap9__48.pngfigure 9.11 : cartes de temps-1 à 4 instants différents


début chapitre

Le pendule paramétrique

Equation du pendule

Le pendule paramétrique non amorti a pour équation générale
MATH

On obtient cette équation lorsque l'un des paramètres du pendule subit une modification au cours du temps, par exemple pour un pendule simple dont le point d'attache oscillerait selon la verticale.

Si h(t) est de la forme MATH on obtient l'équation
MATH

Lorsque MATH cela revient à la balançoire que 2 personnes poussent chacun d'un côté ;il y a résonance paramétrique quand la fréquence de l'excitateur est le double de la fréquence propre du pendule.

Un autre exemple de pendule paramétrique est fourni par une aiguille aimantée placée au centre d'une paire de bobines de Helmhotz dans laquelle on fait passer la superposition d'un courant continu et d'un courant alternatif. L'aiguille est soumise à un couple de forces dont le moment est MATH ;son équation s'écrit, en tenant compte d'un amortissement fluide
MATH

$\theta $ est l'angle de l'aiguille avec la direction du champ fixe d'amplitude $B_{0},$ $J$ le moment d'inertie de l'aiguille, m le moment magnétique de l'aiguille. Cette équation est bien celle d'un pendule paramétrique.

Evolution de la trajectoire

Traçons $\theta $ en fonction de t pour $\omega _{0}=1,$ $\omega =2$ et $\rho =0,2$ dans l'équation du pendule paramétrique :

> restart: with(DEtools):

> phaseportrait([D(x)(t)=u,D(u)(t)=-(1+0.2*sin(2*t))*sin(x)],

[x,u],t=0..200,[[x(0)=0.4,u(0)=0]],stepsize=0.1,

linecolor=black,scene=[t,x]);


chap9__62.pngfigure 9.12 : le mécanisme de l'amplification paramétrique

Lorsque l'amplitude est faible, l'excitateur donne de l'énergie au pendule et l'amplitude de celui-ci augmente. Mais la période du pendule augmente elle aussi avec l'amplitude, ce qui fait que l'excitateur se décale peu à peu et que l'amplitude décroit. On a donc des battements.

En présence d'amortissement, les battements diminuent rapidement et l'amplitude du régime forcé devient constante :

> phaseportrait([D(x)(t)=u,

D(u)(t)=-0.02*u-(1+0.2*sin(2*t))*sin(x)],[x,u],t=0..300,

[[x(0)=0.4,u(0)=0]],stepsize=0.1,linecolor=black,scene=[t,x]);


chap9__63.pngfigure 9.13 : pendule paramétrique avec frottement

La période des oscillations est le double de la période de l'excitateur.

Si l'on augmente $\rho ,$ l'amplitude du pendule augmente, et, du fait des transitoires, l'angle va bientôt dépasser $\pi $ ;le pendule décroche. Pour $\rho =1,1,$ on trouve :

> phaseportrait([D(x)(t)=u,

D(u)(t)=-0.1*u-(1+1.1*sin(2*t))*sin(x)],[x,u],t=42..60,

[[x(0)=0.1,u(0)=0],[x(0)=-0.1,u(0)=0],[x(0)=1,u(0)=0]],

stepsize=0.1,linecolor=black,scene=[t,x]);


chap9__67.pngfigure 9.14 : multistabilité du pendule paramétrique

Pour une même valeur des paramètres, on a trois régimes permanents différents suivant les conditions initiales :des oscillations régulières à une période double de celle de l'excitateur, un régime où le pendule tourne toujours dans le même sens à la fréquence de l'excitateur et un régime identique en sens contraire. Il y a multistabilité.

Si $\rho $ augmente encore, le pendule n'arrive pas à se fixer sur l'un des trois régimes mais saute constamment de l'un à l'autre du fait des transitoires :la trajectoire devient chaotique.

Tracé de l'attracteur

On choisit $\rho =3$ .Le régime est alors chaotique :

> restart: with(DEtools):with(plots):

>phaseportrait([D(x)(t)=u,

D(u)(t)=-0.1*u-(1+3*sin(2*t))*sin(x)],[x,u],

t=0..200,[[x(0)=0.1,u(0)=0]],stepsize=0.1,scene=[t,x]);


chap9__70.pngfigure 9.15 : trajectoire chaotique de l'oscillateur paramétrique

Pour tracer l'attracteur, on change t en $2\pi t$ et x en $2\pi x$ de façon à avoir une période de 1 pour x et de 0,5 pour t. Il faut alors multiplier le second membre par $2\pi =6,28$ :

> p:=phaseportrait([D(x)(t)=u,

D(u)(t)=6.28*(-0.1*u-(1+3*sin(2*t*2*Pi))*sin(x*2*Pi))],[x,u],

t=0..750,[[x(0)=0.1,u(0)=0]],stepsize=0.05,

linecolor=black,scene=[x,u]):

> a:=op(1,op(1,p)):nops(a);

> liste1:=op([]):for i to (nops(a)-1)/10 do

if op(1,op(10*i,a))>0 then liste1:=liste1,

[frac(op(1,op(10*i,a))),op(2,op(10*i,a))] else liste1:=liste1,

[1+frac(op(1,op(10*i,a))),op(2,op(10*i,a))] fi od:

> plot([liste1],style=point);


chap9__74.pngfigure 9.16 : attracteur étrange de la carte de temps-1 du pendule paramétrique

Les ondulations sont là encore le résultat du processus d'étirement-repliement des équations différentielles.

début chapitre

Aiguille dans deux champs magnétiques

On peut également obtenir un mouvement chaotique avec une aiguille aimantée soumise à deux champs magnétiques, l'un fixe et l'autre tournant avec une pulsation $\omega .$ A un instant considéré, l'aiguille fait un angle $\theta $ avec le champ fixe, et, l'angle du champ tournant avec le champ fixe étant $\omega t,$ un angle $\omega t-\theta $ avec le champ tournant ;d'où, le moment de la force magnétique étant MATH on a l'équation, obtenue en projetant sur l'axe perpendiculaire au plan du mouvement
MATH

pour le mouvement non amorti.

L'équation peut s'écrire de forme simplifiée, en ajoutant un frottement fluide, en notant x l'angle à la place de $\theta $
MATH

On prend c = 0,1, $\omega =0,5,$ on multiplie t et x par 2$\pi $ pour avoir une période de 1 pour x et 2 pour t ;$\overset{..}{x}$ est alors divisé par 2$\pi $ et on peut écrire
MATH

La dynamique est complexe. On observe un attracteur chaotique lorsque $\rho $ atteint une valeur proche de 0,7 ;cet attracteur est atteint par une cascade de doublement de période comme on peut le voir en traçant 3 trajectoires :

> restart:with(DEtools):with(plots):

> for i in [0.64,0.67,0.68] do p[i]:=phaseportrait([D(x)(t)=u,

D(u)(t)=2*Pi*(-0.1*u-sin(2*Pi*x)-i*sin(Pi*t-2*Pi*x))],[x,u],

t=80..100,[[x(0)=-0.1,u(0)=0]],stepsize=0.1,scene=[t,x]) od:

> A:=array(1..3,1..1);

> A[1,1]:=p[0.64]:A[2,1]:=p[0.67]:A[3,1]:=p[0.68]:

> display(A);


chap9__89.pngfigure 9.17 : scénario de doublement de période

Pour $\rho =0,4$ (dessin supérieur), la période est de 2, qui est la période du champ tournant.

Pour $\rho =0,67,$ (dessin du milieu), la période est de 4 (de 85 à 89 par exemple) ;c'est le double de la période du champ tournant.

Pour $\rho =0,68,$ la période a encore doublé (de 85 à 94 par exemple) ;elle est de 4 fois la période du champ tournant.

Pour $\rho =0,7,$ le système est chaotique ;il oscille irrégulièrement et saute de temps en temps un ou plusieurs tours ;son attracteur étrange sera l'attracteur A ;traçons la trajectoire, puis l'attracteur :

> phaseportrait([D(x)(t)=u,

D(u)(t)=2*Pi*(-0.1*u-sin(2*Pi*x)-0.7*sin(Pi*t-2*Pi*x))],[x,u],

t=0..200,[[x(0)=-0.1,u(0)=0]],stepsize=0.1,

scene=[t,x],linecolor=red);


chap9__95.pngfigure 9.18 : trajectoire chaotique avec $\rho =0,7$

> p:=phaseportrait([D(x)(t)=u,

D(u)(t)=2*Pi*(-0.1*u-sin(2*Pi*x)-0.7*sin(Pi*t-2*Pi*x))],[x,u],

t=0..1500,[[x(0)=-0.1,u(0)=0]],stepsize=0.1,scene=[x,u]):

> a:=op(1,op(1,p)):nops(a);

> liste1:=op([]):for i to (nops(a)-1)/20 do

liste1:=liste1,[frac(op(1,op(20*i,a))),op(2,op(20*i,a))] od:

> liste2:=op([]):for i to nops([liste1]) do

if op(1,op(i,[liste1]))>0.7 then liste2:=liste2,

[op(1,op(i,[liste1]))-1,op(2,op(i,[liste1]))]

else liste2:=liste2,[liste1][i] fi od:

> plot([liste2],style=point);


chap9__97.pngfigure 9.19 : attracteur de la boussole avec 2 champs magnétiques pour $\rho =0,7$

Pour une valeur de $\rho $ voisine de 0,725, la boussole a tendance à s'accrocher sur le champ tournant ;l'attracteur grandit brusquement et devient l'attracteur B :il subit une crise :


chap9__100.pngfigure 9.20 : attracteur pour $\rho =0,725$

Lorsque $\rho $ continue à croître, on a une succession de trajectoires régulières et chaotiques avec un attracteur de type B.

début chapitre

Trajectoires d'une boule de Billard

La trajectoire d'une boule qui rebondit entre deux parois obéit à la loi de Descartes de la réflexion à chaque choc (angle de réflexion égal à l'angle d'incidence). Cette loi simple peut aboutir à une trajectoire chaotique si la boule rebondit entre un carré et un cercle de même centre à l'intérieur du carré.

Billard entre deux cercles

Commençons par le cas simple d'une boule qui rebondit entre deux cercles de même centre.

Les 2 cercles ont leur centres en O et sont de rayon 1 et 2. Le point de départ est le point $M_{0}$ $(x_{0}=2,$ $y_{0}=0)$ et le premier point $M_{1}$ est sur le cercle de rayon 1 et sur la droite inclinée de $\alpha =$13$^{\text{o}}$ sur l'horizontale, de coordonnées (MATH MATH Le point suivant $M_{2}$ $(x_{2},$ $y_{2})$ sur le grand cercle sera symétrique de $M_{0}$ par rapport au rayon passant par $M_{1}.$ On a donc MATH k étant une constante,soit
MATH

De plus la distance $OM_{2}=2.$ Pour éliminer le symétrique de $M_{0} $ par rapport au centre qui satisferait les équations, on ajoute la condition $M_{0}M_{2}<1.$ Les points pairs sont obtenus de cette façon, les points impairs ont une distance au centre égale à 1.

> restart:with(plots):with(plottools):

> alpha:=evalf(13*Pi/180);x0:=2;y0:=0;

x1:=cos(alpha);y1:=sin(alpha);(point initial)

> liste:=line([x0,y0],[x1,y1]): for i to 22 do if frac(i/2)<>0 then (si n impair)

x2:='x2':y2:='y2':sol:=solve({(x2+x0)/x1=(y2+y0)/y1,

x2^2+y2^2=4,(x2-x0)^2+(y2-y0)^2<1},{x2,y2});

assign(sol): liste:=liste,line([x1,y1],[x2,y2]):x0:=x1:y0:=y1:

x1:=x2:y1:=y2

elif frac(i/2)=0 then (si n pair)

x2:='x2':y2:='y2':sol:=solve({(x2+x0)/x1=(y2+y0)/y1,

x2^2+y2^2=1,(x2-x0)^2+(y2-y0)^2<1},{x2,y2});

assign(sol): liste:=liste,line([x1,y1],[x2,y2]):

x0:=x1:y0:=y1:x1:=x2:y1:=y2 fi od:

> liste1:=liste:

> liste2:=liste:(on recommence avec MATH)

> p1:=circle([0,0],1,thickness=3):p2:=circle([0,0],2,thickness=3):

> display(p1,p2,liste1,liste2,scaling=constrained,color=red);


chap9__121.pngfigure 9.21 : billard entre 2 cercles

L'erreur est multipliée par 2 à chaque réflexion, mais le mouvement n'est pas chaotique, car la corrélation entre les 2 mouvements voisins à l'origine subsiste. Recommençons avec un billard carré.

Billard carré

Cette fois, la boule rebondit, en suivant la loi de réflexion de Descartes, entre un cercle de rayon 1 et un carré de même centre de côté 2.

La boule part du point initial $(x_{0},y_{0})$ ;elle vient taper sur le cercle si la valeur absolue de la distance de l'origine à la droite initiale est inférieure au rayon du cercle ;cette distance d est donnée par $d=p\cos (\alpha )$, où p est l'ordonnée à l'origine de la droite d'équation MATH soit MATH

Si donc la boule vient taper le cercle, le point d'impact est $(x_{1},y_{1})$ commun à la droite et au cercle tel que la valeur absolue de la différence des abcisses $x_{1}-x_{0}$ soit minimale, pour écarter l'autre intersection.

La trajectoire réfléchie sera alors symétrique de la trajectoire incidente par rapport à la normale, c'est-à-dire le rayon passant par $(x_{1},y_{1}).$ Si $\beta $ est l'angle de ce rayon avec l'horizontale MATH l'angle de la trajectoire réfléchie avec l'horizontale sera MATH $(x_{2},y_{2})$ sera l'intersection de la droite réfléchie passant par $(x_{1},y_{1})$ de pente $\tan (\delta )$ avec le premier côté du carré rencontré.

Pour cela, on cherche les 4 intersections de la droite avec les 4 côtés ou prolongements des côtés ;l'intersection correcte sera celle qui coupe le côté et non son prolongement et dans le sens du rayon, c'est-à-dire que la projection du rayon réfléchi sur le rayon doit être > 0.

Reste le cas où la trajectoire ne coupe pas le cercle. Alors l'intersection se fait sur le carré et non sur son prolongement au point différent du point initial. L'angle réfléchi est alors l'opposé de l'angle d'incidence.


chap9__135.pngfigure 9.22

Tout cela donne le programme suivant :

> restart:with(plottools):with(plots):

> r:=1:a:=2:x0:=2:y0:=0:alpha:=evalf(-15*Pi/180);

> liste:=op([]): for i to 100 do: p:=y0-tan(alpha)*x0:

d:=abs(p*cos(alpha)):(premier cas, la particule tape le cercle)

if d<r then sol:=solve({y=tan(alpha)*x+p,x^2+y^2=r^2},{x,y}):

c:=subs(sol[1],x): e:=subs(sol[2],x): if abs(x0-c)<abs(x0-e)

then x1:=c:y1:=subs(sol[1],y) else x1:=e:y1:=subs(sol[2],y)

fi:(on écarte l'autre intersection)

beta:=arctan(y1/x1):delta:=-alpha+2*beta:

eq:=y=tan(delta)*(x-x1)+y1:c1:=subs(x=a,rhs(eq)):

c2:=subs(x=-a,rhs(eq)):d1:=solve(subs(y=a,eq),x):

d2:=solve(subs(y=-a,eq),x):(l'intersection correcte est celle qui coupe le côté dans le sens du rayon)

if abs(c1)<a and (a-x1)*x1+(c1-y1)*y1>0 then x2:=a:y2:=c1

elif abs(c2)<a and (-a-x1)*x1+(c2-y1)*y1>0 then x2:=-a:

y2:=c2:elif abs(d1)<a and (d1-x1)*x1+(a-y1)*y1>0

then x2:=d1:y2:=a

else x2:=d2:y2:=-a fi:liste:=liste,line([x0,y0],[x1,y1]),

line([x1,y1],[x2,y2]): x0:=x2:y0:=y2:alpha:=-delta:

else

eq:=y=tan(alpha)*x+p:c1:=subs(x=a,rhs(eq)):

c2:=subs(x=-a,rhs(eq)):d1:=solve(subs(y=a,eq),x):

d2:=solve(subs(y=-a,eq),x):

if abs(c1)<a and x0<>a then x2:=a:y2:=c1

elif abs(c2)<a and x0<>-a then x2:=-a:y2:=c2

elif abs(d1)<a and y0<>a then x2:=d1:y2:=a

else x2:=d2:y2:=-a fi: liste:=liste,line([x0,y0],[x2,y2]):

x0:=x2:y0:=y2:alpha:=-alpha fi od :

> p1:=circle([0,0],1):p2:=rectangle([-2,2],[2,-2]):

> display(p1,p2,liste,scaling=constrained);


chap9__136.pngfigure 9.23 : trajectoire chaotique du billard carré

La trajectoire de la boule est chaotique et soumise à la sensibilité aux conditions initiales ;si l'on fait partir 2 boules du même point initial avec des angles très voisins, les trajectoires sont très rapidement déconnectées.

début chapitre

retour index