TangenteX.com

Systèmes non-linéaires et chaos

Un système non-linéaire, encore appelé "système dynamique", modélise un phénomène à l'aide d'équations non linéaires. Il peut s'agir d'équations différentielles, aux dérivées partielles ou encore des équations aux différences (ce sont des équations discrètes comme la définition de suites récurrentes). Une équation non linéaire contient dans un de ses membres une composante non linéaire, comme par exemple un terme quadratique ou trigonométrique. Exemple connu sur TangenteX d'équation non linéaire : \( \dfrac{d^2 \theta}{dt^2} + \omega_0^2 \sin(\theta) = 0 \), qui contient un sinus, terme non linéaire, car je rappelle que sin(a + b) n'est pas égal à sin(a) + sin(b) !

Les systèmes non linéaires sont extrémement courants dans tous les domaines des sciences physiques et même en sciences humaines comme l'économie ou la sociologie. Les phénomènes qui répondent strictement à des systèmes linéaires sont très exceptionnels. Et pourtant, les physiciens ont une forte tendance à essayer de linéariser les systèmes. D'ailleurs en prépa et en licence, la plupart de systèmes étudiés sont linéarisés : pourquoi ?

Les systèmes linéaires, et les équations les modélisant sont d'étude simple, du fait de la propriété de linéarité : on peut décomposer le problème, trouver des solutions particulières puis déterminer la solution générale en combinant linéairement des solutions particulières. Impossible lorsqu'il s'agit d'équations non linéaires.

Cette particularité a longtemps freiné l'étude des systèmes non linéaires. Les mathématiciens qui les abordèrent le faisaient comme des naturalistes : ils répertoriaient les différentes formes de systèmes et leurs comportements. Le plus illustre de ces mathématiciens est sans doute Henri Poincaré. Tout changea dans la seconde moitié du XXième siècle, avec l'apparition des ordinateurs. Il devint possible de résoudre numériquement les systèmes non linéaires, qui ne possèdent pas, pour la plupart, de solution analytique. Et l'on vit apparaître les premiers résultats à propos des systèmes dynamiques et du chaos, par exemple avec Lorenz en 1963.

Pourquoi associer le chaos aux systèmes non linéaires ? Comment caractériser le chaos ?

La notion de chaos fut théorisée par Poincaré dans son étude sur la stabilité des systèmes non linéaires. Le comportement chaotique est intimement lié au caractère non linéaire d'un système, il est impossible de produire du chaos dans un système linéaire. Le comportement chaotique se retrouve dans une grande diversité de phénomènes physiques, de la turbulence dans les fluides à l'arythmie cardiaque, des interactions proies-prédateurs aux cavités laser.

Mais attention à ne pas confondre système chaotique et système stochastique (aléatoire). Le chaos n'est pas aléatoire. Certes, vu de loin, on pourrait les confondre. La marche aléatoire d'un ivrogne peut sembler chaotique, et le langage courant supporte souvent la confusion. Il existe cependant une différence essentielle entre le stochastique et le chaotique : un système chaotique est déterministe, alors qu'un système stochastique ne l'est pas. L'information nécessaire pour décrire entièrement un système stochastique est infinie. L'information pour décrire entièrement un système chaotique tient dans une équation ! Cette différence est fondamentale.

Toutefois, les systèmes stochastiques peuvent donner lieu à l'émergence de comportements déterministes, en particulier lors de changement d'échelle. A l'échelle microscopique, tous les systèmes sont gouvernés par les lois de probabilité quantiques. Pourtant, à l'échelle macroscopique, il existe des systèmes déterministes ! A ce propos, je vous signale un livre de David Ruelle "Hasard et chaos" ([1] voir bibliographie), dont j'encourage vivement la lecture.

Une des caractéristiques déterminantes d'un système chaotique est sa sensibilité aux conditions initiales (SCI), qui le rend non prédictible en pratique. Son caractère déterministe fait que l'on devrait toujours savoir prédire le comportement du système en étudiant ses équations et ses conditions initiales (il existe des théorèmes en mathématique qui nous le garantissent). Sauf que dans la pratique, nous sommes incapables de définir les conditions initiales avec suffisamment de précision pour calculer l'évolution du système sans se tromper au délà d'une durée plus ou moins grande.

Un système chaotique est déterministe et non prédictible. Un système stochastique est non prédictible, aussi, mais n'est pas déterministe.

Notons qu'il n'existe pas de définition précise de ce qu'est le chaos en dynamique des systèmes, mais que l'on sait exhiber différentes propriétés qui le caractérisent : l'aspect de son portait de phases et de sa section de Poincaré, son spectre de puissance ou en encore l'existence d'un exposant de Lyapunov positif. Nous verrons cela avec les simulations.

Le modèle du pendule physique

Pour aller à la découverte du chaos, j'ai choisi un modèle bien connu de tous : le pendule. Le modèle du pendule s'organise en poupées russes : on part du plus simple, le pendule simple linéaire, pour évoluer vers le pendule amorti puis le pendule forcé, et enfin le modèle non-linéaire. Ce cheminement se fait en complexifiant le modèle pour qu'il colle de plus en plus à la réalité. J'ai traité ces différents aspects dans la page de TangenteX consacrée aux oscillateurs.

Le pendule physique

La littérature américaine nomme "pendule physique" le modèle du pendule non linéaire, amorti et forcé. Parce que cette notation est simple et assez parlante, je la conserverai. Nous avons maintes fois rencontré dans les pages de TangenteX l'équation différentielle qui modélise le pendule simple non linéaire :

\( \dfrac{d^2 \theta}{dt^2} + \omega_0^2 \sin(\theta) = 0 \)

La fonction \( \theta (t) \) décrivant la variation angulaire du pendule avec la verticale. La pulsation propre du pendule est notée \( \omega_0 \). Sous cette forme, le système ne dissipe pas d'énergie, on dit qu'il est "hamiltonien" ou encore "conservatif".

Ajoutons un peu de réalisme à ce modèle en lui ajoutant une force de trainée aérodynamique qui va consommer de l'énergie du pendule, ce qui va avoir pour effet de l'amortir. Dans ce modèle, je suppose que cette force est proportionnelle à sa vitesse angulaire (l'amortissemnt est faible), ce qui me donne l'équation :

\( \dfrac{d^2 \theta}{dt^2} + \beta  \dfrac{d \theta}{dt} + \omega_0^2 \sin(\theta) = 0 \)

Le coefficient \( \beta \) est le coefficient d'amortissement. Cette équation définit un système non-linéaire, par la présence du sinus, et dissipatif, par la présence du coefficient d'amortissement \( \beta \).

Enfin, je vais exciter mon pendule avec un signal périodique, afin de lui restituer l'énergie qu'il perd par frottement. J'ai choisi un signal sinusoïdal, d'amplitude \( A_F \) et de pulsation \( \omega_F \), ce qui me donne l'équation :

\( \dfrac{d^2 \theta}{dt^2} + \beta  \dfrac{d \theta}{dt} + \omega_0^2 \sin(\theta) = A_F \sin(\omega_F t) \)            (1)

Par un changement de variables, j'adimensionne l'équation (1) afin de rendre son étude plus facile, ce qui me donne :

\( \dfrac{d^2 \theta}{dt^2} + \gamma  \dfrac{d \theta}{dt} + \sin(\theta) = A \sin(\omega_D t) \)            (2)

L'équation (2) définit un système non linéaire de dimension 2, non autonome (il dépend du temps), dissipatif et forçé. L'équation (2) s'exprime donc sous la forme du système différentiel :

\( \left\{ \begin{array}{l r} \dfrac{d \theta}{dt} = y \\ \dfrac{d y}{dt} = -\gamma y - \sin(\theta) + A\sin(\omega_D t) \end{array} \right. \)

Et là, je vais faire appel à deux théorèmes fondamentaux des systèmes différentiels :

  • Le théorème de Cauchy-Lipschitz, qui me garantit que pour un jeu de conditions initiales données, 2 ici car mon équation différentielle est d'ordre 2, il n'existe qu'une seule solution à l'équation (2).
  • Le théorème de Poincaré-Bendixson, qui me dit que pour un système différentiel d'ordre inférieur à 3, il ne peut y avoir d'apparition du chaos. Dit avec les mains, dans un plan, les trajectoires n'ont pas la place pour s'emmêler sans se couper.

Le premier théorème garantit le déterminisme du système et m'arrange bien. Le second est beaucoup plus ennuyeux dans la mesure où mon objectif est de faire apparaître le chaos dans mon système, et que ce théorème dit que c'est impossible. Cependant, on constate le chaos expérimentalement : pourquoi ? Parce que sous ses apparences, (2) peut être écrite sous la forme d'un système de dimension 3. Procédons à un changement de variable, en posant \( z = \omega_D t \) et ré-écrivons le système différentiel :

\( \left\{ \begin{array}{l r} \dfrac{d \theta}{dt} = y \\ \dfrac{d y}{dt} = -\gamma y - \sin(\theta) + A\sin(z) \\  \dfrac{dz}{dt} = \omega_D \end{array} \right. \)

J'ai introduit la variable z, qui est la phase du terme de forçage dans l'équation (2). Nous voilà donc en présence d'un système à trois variables \( \theta, \dfrac{d \theta}{dt}, z \) qui permet de modéliser l'apparition du chaos. Il possède également trois paramètres \( \gamma, A, \omega_D \), dont chacun peut constituer un paramètre de contrôle du système.

Dans la suite de nos manipulations, nous allons faire varier l'amplitude A de forçage du pendule pour étudier son comportement. A sera le paramètre de contrôle de notre système, sa variation provoquant des changements qualitatifs dans l'évolution du système.

Son modèle numérique

Le modèle numérique de notre système est codé dans la fonction :

def SYS(Theta,t):
    thetaPoint = Theta[1]
    dthetaPoint = -gamma*Theta[1] - sin(Theta[0]) + A*sin(wD*t)
    return array([thetaPoint, dthetaPoint])

Theta est une liste de bipoints, qui contient pour chaque pas de temps dans la première colonne l'angle du pendule avec la verticale et dans la deuxième colonne la vitesse angulaire.

Dans cette implémentation, les paramètres gamma, A et wD sont déclarés en variables globales des scripts. J'ai fait une exception dans le script Bifurcation.py, dans lequel je passe A en paramètre d'appel de la fonction SYS.

Vous remarquerez dans les différents scripts la fonction Recalage. Cette fonction permet de ramener la valeur de l'angle de l'oscillation dans le domaine de variation \( [-\pi, \pi] \), car il ne vous aura pas échappé que notre pendule pouvait faire une ou plusieurs rotations autour de son axe. Son code est simple et illustre le mécanisme :

def recalage(theta):
    for i, angle in enumerate(theta[:,0]):
        while angle > pi:
            angle -= 2.0*pi
        while angle < -pi:
            angle += 2.0*pi
        theta[i,0] = angle   
    return theta

Les paramètres de base du modèle

Vous pouvez bien entendu modifier tous les paramètres du modèle, la pulsation propre du pendule, le coefficient d'amortissement, les paramètres de forçage. Mais aussi le nombre de périodes affichées ou la fréquence d'échantillonnage.

Cependant, pour que l'étude soit plus simple, j'ai fixé certains paramètres aux valeurs suivantes :

gamma = 0.2          # coefficient d'amortissement
Fe = 1000.0          # fréquence d'échantillonnage du signal
NbPer = 10           # nombre de périodes de forçage affichées

J'ai aussi fixé les conditions initiales du modèle pour la première partie de la simulation :

theta0 = pi/10.      # angle initial
thetaPoint0 = 0.0    # vitesse angulaire initiale

Les paramètres que je ferai varier sont AF et wF, respectivement l'amplitude et la pulsation du signal de forçage.

Enfin, deux paramètres vous permettent de déclencher ou non le calcul du spectre et de la section de Poincaré dans le script PenduleChaos.py, car ces calculs peuvent être longs, surtout celui de la section de Poincaré :

Aff_Spectre = True
Aff_Poincare = True

Si vous ne voulez pas déclencher un de ces deux calculs, ou les deux, fixez la valeur du paramètre concerné à False.

Les outils d'analyse

Il existe plusieurs outils mathématiques pour mettre en évidence et analyser le comportement d'un système dynamique. J'en expose quelques uns ici, que j'ai considérablement simplifié. En particulier, je n'aborderai pas les mathématiques sous-jacentes. J'indiquerai à la fin de cette page quelques ouvrages qui intéresseront ceux qui souhaitent creuser le sujet, en particulier l'ouvrage de Robert Roussarie ([4] voir bibliographie).

Le portait de phase

Vous avez sans doute l'habitude de représenter l'évolution d'un système mécanique comme notre pendule physique par une courbe d'évolution, c'est à dire la variation d'une variable représentative du système, l'angle \( \theta \) généralement pour notre pendule, en fonction du temps. Il existe d'autres représentations possibles, dont la plus connue est le portait de phase.

On appelle "trajectoire", la courbe \( C(x, \dfrac{dx}{dt}) \) pour des conditions initiales données. Dans le cas de notre pendule, les conditions initiales sont \(( \theta_0, \dfrac{d \theta}{dt}_{t = 0} \), l'angle et la vitesse angulaire initiales. Chaque trajectoire est donnée par la courbe \( C(\theta, \dfrac{d \theta}{dt})\) pour un couple de données initiales donné.

Le portait de phase d'un système est un ensemble de trajectoires du système, exprimées à l'aide des quantités conjuguées \( \theta \) et \( \dfrac{d \theta}{dt}\)  On le trace en résolvant le système différentiel pour un ensemble de conditions initiales données.

Sur un portait de phase, nous pouvons distinguer les attracteurs et les bassins d'attraction des trajectoires. Ces noms sont intuitifs : ils désignent des points et des régions du plan de phase vers lesquelles convergent les trajectoires lorsque t tend vers l'infini. L'existence de ces attracteurs permet de juger intuitivement de la stabilité ou de l'instabilité du système. Nous verrons lors des simulations ce que nous pouvons en tirer.

La section de Poincaré

Le portrait de phase d'un système dynamique peut contenir des trajectoires très compliquées et embrouillées, surtout lorsqu'il présente un comportement chaotique. Poincaré l'a constaté lors de son étude de la dynamique du système solaire. Aussi, a-t-il inventé, en 1892, un outil simple et universel : la section de Poincaré ou coupe de Poincaré (les anglophones utilisent "Poincaré plot"). N'oublions pas qu'il faisait ses calculs à la main...

Nous avons vu que notre système dynamique est de dimension 3, c'est à dire que ses trajectoires évoluent dans l'espace. Pour étudier cette évolution, Poincaré a imaginé de couper en tranches cet espace par une surface, la section de Poincaré,  et d'examiner la forme de l'intersection de cet espace de dimension 3 avec une surface de dimension 2, c'est à dire un ensemble de points de la surface. D'ailleurs, le principe est généralisable à des variétés de dimensions supérieures. On réduit donc la dimension du système à étudier, en transformant l'infinité de points d'une trajectoire de phase en un point sur la section de Poincaré.

Pratiquement, on utilise ce que l'on appelle la méthode stroboscopique. On prend une "photo", c'est à dire une coupe par une surface de l'ensemble des trajectoires du système, à intervalle temporel régulier. Il reste à fixer la valeur de cet intervalle. Comme lorsqu'on réalise l'étude d'un mouvement en l'éclairant avec une lumière stroboscopique, je vais fixer cette valeur égale à la période du signal de forçage \( T_F = \dfrac{2 \pi}{\omega_D}\).

Nous verrons plus loin les informations que nous pouvons tirer de la section de Poincaré de notre système "Pendule physique", pour différentes valeurs de la période de forçage du système.

Du point de vue programmation, le calcul de la section de Poincaré d'un système est simple, mais assez consommatrice de temps CPU. Voici le code que j'utilise :

def SectionPoincare(sys,c0,perF,n):
    t = linspace(0, n*perF,n**2)
    SP = odeint(sys,c0,t)
    SP = recalage(SP)
    x = [SP[n*i,0] for i in range(n)]
    y = [SP[n*i,1] for i in range(n)]
    return x,y

Le principe consiste à calculer les trajectoires du système avec odeint() puis de recaler les variations angulaires dans l'intervalle \( [- \pi, \pi] \). Ceci fait, je calcule les coordonnées x et y des points d'intersection des trajectoires stockées dans le tableau SP, pour n surfaces coupant les trajectoire tous les perF (la période de forçage) intervalles de temps.

Le spectre de puissance

Le spectre de puissance d'un système chaotique comporte une partie continue et une bande large. Il est assez typique, en faisant apparaître une forêt de pics plus ou moins importants, alors que le spectre de puissance d'un système périodique ou quasi-périodique comporte une seule raie (un pic) avec quelques raies secondaires, correspondant aux harmoniques. L'aspect du spectre de puissance d'un système est donc un indice important d'apparition du chaos.

Le spectre de puissance ne permet cependant pas, en général, de distinguer un système chaotique déterministe d'un système stochastique générant un bruit "blanc" (un bruit dont la densité spectrale de puissance est la même en moyenne pour toutes les fréquences). Toutefois, dans les systèmes chaotiques forçés comme notre pendule, nous pouvons repérer des pics correspondant à la fréquence de forçage et à ses harmoniques, qui ne peuvent pas exister, par définition, dans le bruit blanc généré par un système stochastique.

Pour étudier un système chaotique, il faut associer l'analyse du spectre de puissance du système à celle de sa section de Poincaré.

Le code nécessaire au calcul du spectre de puissance est bien connu, nous l'avons abordé dans cette page de TangenteX :

def Spectre(theta,Te):
    N = theta[:,0].size
    theta_FFT = abs(fft(theta[:,0]))/N
    theta_freq = fftfreq(N,Te)
    theta_FFT = theta_FFT[0:len(theta_FFT)//2]
    theta_freq = theta_freq[0:len(theta_freq)//2]
    return theta_freq,theta_FFF

L'exposant de Lyapunov

Une des caractéristiques d'un système chaotique est sa sensibilité aux conditions initiales. Cela signifie que pour des conditions initiales très proches, les trajectoires divergent plus ou moins rapidement. L'exposant de Lyapunov (en fait les exposants, mais je me limiterai à l'exposant majeur) permet une estimation de cette vitesse de divergence.

Pour mettre en évidence l'exposant de Lyapunov, je vais donc calculer deux trajectoires du système "Pendule physique" pour les mêmes paramètres mais avec des conditions initiales différentes, d'une valeur très faible. Puis je trace la courbe semi-logarithmique de la différence de variable angulaire en fonction du temps, pour constater leur évolution respective.

Aleksandr Lyapunov a établi que la variation de cette différence répond à la loi \( \Delta \theta(t)   \approx C e^{\lambda t}\), avec \( \lambda \) l'exposant de Lyapunov. Si nous traçons la courbe des différences de trajectoires en coordonnées semilog, alors nous pouvons estimer l'exposant de Lyapunov en mesurant la pente de cette courbe.

La valeur de l'exposant de Lyapunov est un des critères de détermination du comportement chaotique du système. De la définition de l'exponentielle, nous en déduisons immédiatement qu'un comportement chaotique, dans son aspect sensibilité aux conditions initiales, implique que \( \lambda \) doit être positif. En effet, pour \( \lambda \gt 0 \), la différence entre les deux trajectoires croît de façon exponentielle.

Le script Lyapunov.py permet le calcul et la visualisation de l'exposant de Lyapunov.

Le diagramme de bifurcation

Un autre des critères de mise en évidence du chaos dans un système s'observe sur le tracé de son diagramme de bifurcation.

Une bifurcation est un changement brutal dans la dynamique d'un système, qui survient lors de l'évolution d'un paramètre du système, par exemple l'amplitude de forçage de notre pendule physique. On appelle le point auquel survient ce changement un point critique du système.

Les bifurcations sont des phénomènes typiques des systèmes non linéaires. Dans un système linéaire, les effets sont proportionnels aux causes : petites causes, petits effets, grandes causes, grands effets. Ce n'est plus vrai pour les systèmes non linéaires, qui présentent des variations qui peuvent être sans rapport avec les causes initiales. C'est le fameux effet "papillon".

Il existe une théorie mathématique assez vaste des bifurcations, dont vous trouverez une introduction dans [4] tome 2, par exemple. Ici, je me contenterai de mettre en évidence visuellement l'existence de ces bifurcations.

Pour cela, j'exploiterai une caractéristique du comportement non linéaire. A l'approche d'une bifurcation la distance sur le plan de phase des trajectoires devient anarchique. Je vais donc calculer dans un domaine de variation donné de l'amplitude de forçage, la trajectoire du système. Après avoir supprimé l'intervalle de temps qui correspond au régime transitoire du système, je calcule pour les n premiers points de la trajectoire, la distance séparant la position angulaire et la vitesse angulaire. J'ai fixé n à 200 avec mes valeurs de paramètres. Puis je trace simplement, pour chaque valeur de l'amplitude de forçage, la répartition des distances. Nous verrons plus loin ce que cela donne.

Le code de calcul utilisé est :

for A in AF:
    SP = odeint(SYS,C0,t,args=(A,))
    mask_transitoire = t > 10*TF
    SP = SP[mask_transitoire]
    r[i] = hypot(SP[200,0],SP[200,1])
    i += 1

Le script Bifurcation.py permet le calcul et la visualisation du diagramme de bifurcation.

Les programmes de simulation

Pour mettre en oeuvre les différents outils que je viens de présenter, j'ai écrit trois scripts Python :

  • PenduleChaos.py, qui trace pour un jeu donné de paramètres, la courbe d'évolution du système, son portait de phase, son spectre de puissance et sa section de Poincaré.
  • Bifurcation.py, qui construit le diagramme de bifurcation du système sur un intervalle donné d'amplitudes de forçage.
  • Lyapunov.py, qui calcule et visualise la distance qui sépare deux courbe d'évolution pour des conditions initiales très proches. Il permet la mise en évidence de l'exposant de Lyapunov.

Ces scripts Python sont téléchargeables depuis la bibliothèque de codes de TangenteX.

Dans la suite de cette page, j'utiliserai ces outils pour analyser le comportement du système "Pendule physique" en faisant varier les différents paramètres du modèle.

Le régime amorti

Dans l'équation (2), nous pouvons annuler le terme de droite, en fixant l'amplitude de forçage A = 0. Ainsi, nous sommes en présence d'un système non linéaire (le terme en sinus est toujours présent) dissipatif pur. La dissipation d'énergie est proportionnelle au coefficient d'amortissement \( \gamma \) Nous avons déjà rencontré ce système dans la page de TangenteX consacrée aux oscillateurs.

Pour avoir un point de comparaison avec les simulations suivantes, je vais tracer la courbe d'évolution, le portait de phase, le spectre de puissance et la section de Poincaré du pendule physique dissipatif.

Les paramètres du mode

Voici la valeur des paramètres du pendule physique choisis pour la simulation :

wD = 0.8             # pulsation de forçage
A  = 0.0             # amplitude de forçage

Pulsation de forçage : 0.800000
Fréquence de forçage (Hz) : 0.127324
Amplitude de forçage : 0.000000

Evolution et portrait de phase

L'exécution du script PenduleChaos.py restitue l'évolution et portrait de phase suivants :

Trajectoire et portrait de phase sans forçage

Sur la courbe d'évolution, ou "time series" comme disent les anglophones, nous observons une décroissance des oscillations du pendule, décroissance plus ou moins rapide selon la valeur du coefficient d'amortissement.

L'énergie initiale du pendule est amenée par l'énergie potentielle dûe à l'angle initial non nul. Lorsqu'on lâche le pendule, cette énergie décroît pour revenir à la valeur de l'énergie mécanique du pendule au repos, c'est à dire pour un angle avec la verticale nul. C'est ce que l'on observe sur la courbe d'évolution.

Le portait de phase illustre le même phénomène sous un autre point de vue. La trajectoire de phase part du point initial pour lequel l'angle est égal à l'angle initial et la vitesse initiale nulle, c'est à dire nos conditions initiales. Puis la trajectoire parcoure le plan de phase en spiralant vers le point où l'énergie du pendule est la plus faible, de coordonnées (0,0) sur le plan de phase. La trajectoire semble "tomber" dans ce point, qui est le point de moindre énergie, en fait le fond du puit de potentiel du pendule. On appelle ce genre de point, comme notre point (0,0) un attracteur du système.

Spectre de puissance et section de Poincaré

Voici le spectre de puissance et la section de Poincaré pour un régime dissipatif pur :

Spectre et SP sans forçage

Sur le spectre de puissance, nous observons la raie unique qui correspond à la fréquence propre du pendule physique. Sa pulsation propre \( \omega_0 \) est égale à 1.0 et donc sa fréquence propre est égale à environ 0.16 Hz.

Sur la section de Poincaré, nous observons le point initial de la trajectoire, en haut à droite, puis quelques points intermédiaires. Les points se concentrent autour de la position d'équilibre du pendule, à angle nul et vitesse angulaire nulle. C'est cohérent avec le plan de phase.

Le régime forcé périodique

Ajoutons maintenant à notre pendule amorti une source d'énergie extérieure, sous la forme d'une excitation par un signal sinusoïdal d'amplitude et de pulsation paramétrables. On appelle ce type de signal un signal de forçage. Dans la suite de mes manips, je ne fais varier que l'amplitude du signal de forçage, mais bien sur, vous pourrez faire aussi varier sa pulsation.

Les paramètres du mode

Je vais pour le moment choisir une faible amplitude de forçage, en utilisant les paramètres suivants :

wD = 0.8             # pulsation de forçage
A   = 0.5            # amplitude de forçage

Pulsation de forçage : 0.800000
Fréquence de forçage (Hz) : 0.127324
Amplitude de forçage : 0.500000

Vous noterez que la fréquence de forçage, 0.13 Hz, est proche de la fréquence propre du pendule physique, environ 0.16 Hz. Le rapport entre ces deux fréquences est aussi un élément déclencheur du chaos. Je ne le traite pas dans cette page, mais vous pourrez l'expérimenter en fixant l'amplitude de forçage et en faisant varier sa pulsation. J'ai choisi ici de fixer la pulsation de forçage et de faire varier son amplitude.

Evolution et portrait de phase

Pour une valeur faible d'amplitude de forçage, voici la courbe d'évolution et le portait de phase obtenus :

Trajectoire et portrait de phase AF = 0.5

En mode forcé, le pendule reçoit de l'énergie extérieure qui compense la perte d'énergie par amortissement. Sur la courbe d'évolution, nous observons qu'après un régime transitoire qui dure une quarantaine de secondes, soit environ 6 périodes de forçage, le pendule oscille à la période de forçage, avec une amplitude d'oscillation constante.

Le portait de phase du pendule forcé périodique est très différent de celui du pendule amorti. Nous y observons un attracteur, mais qui n'est plus un puit de potentiel mais plutôt un ensemble de trajectoires plus ou moins embrouillées : on parle de "cycle limite". Il est le témoin du régime transitoire puis du régime forcé. On retrouve la forme de la trajectoire du pendule simple non linéaire sans amortissement.

Aucune trajectoire ne part vers l'infini, elles restent toutes dans le domaine de l'attracteur. Physiquement, cela signifie que l'amortissement du pendule permet encore de compenser l'énergie de forçage.

Spectre de puissance et section de Poincaré

Spectre et SP AF = 0.5

Le spectre de puissance du signal montre une raie principale dont la fréquence correspond à la fréquence de forçage du pendule. Mais on observe aussi quelques raies secondaires qui correspondent au régime transitoire en début d'évolution du système.

La section de Poincaré montre un regroupement de points d'intersections qui correspond à l'attracteur constitué par la trajectoire de forçage du pendule. Les poinst isolés correspondent au régime transitoire.

L'apparition du chaos

Nous allons maintenant augmenter l'amplitude de forçage du pendule, paramètre que j'ai choisi comme paramètre de contrôle, en laissant fixes les autres paramètres et les conditions initiales du mouvement.

Diagramme de bifurcation

J'ai indiqué plus haut que le tracé du diagramme de bifurcation d'un système était un outil puissant de détection de l'apparition du chaos. J'ai donc utilisé le script Bifurcation.py pour tracer les bifurcations éventuelles de notre pendule physique en faisant varier l'amplitude de forçage entre 0.0 et 1.5, tous paramètres constants par ailleurs. Voici ce que j'obtiens :

Diagramme de bifurcation

Qu'observe-t-on ? Sur certaines parties du diagramme, par exemple pour les valeurs d'amplitude de forçage comprises entre 0.0 et 0.5, les distances des trajectoires sur le plan de phase sont très petites et stables. Le diagramme fait alors apparaître des plateaux, qui sont le témoin d'un régime périodique. Pour d'autres valeurs de l'amplitude de forçage, vous constatez que les distances sont très variables. Le régime n'est plus périodique mais chaotique. Plus la répartition des points est étendue sur l'axe des ordonnées, plus le régime est chaotique. Il l'est, par exemple, pour A = 1.2, valeur d'amplitude que je choisirai pour la suite des expériences.

A noter sur ce diagramme :

  • La transition entre régime périodique et chaotique est très brusque, à partir d'un point que l'on nomme "point critique". Pour une très petite variation de A, le régime passe de périodique à chaotique. On appelle ce changement de régime un "doublement de période". C'est l'origine du terme "bifurcation", qui emporte l'idée d'un changement brutal.
  • Dans le domaine de variation de A, on observe plusieurs transitions périodique/chaotique. Seul le diagramme de bifurcation permet la visualisation de ce comportement, qui est important dans l'analyse et l'ingénierie du système. Si vous voulez construire un pendule physique forcé en régime périodique, ce qui est généralement le but, vous choisirez une valeur de A pour laquelle le diagramme présente  un plateau.

Evolution et portrait de phase

Pour la valeur A = 1.2 de l'amplitude de forçage, voici la courbe d'évolution et le portait de phase obtenus :

Trajectoire et portrait de phase AF = 1.5

Comme vous le constatez sur l'évolution du système, son régime n'a plus rien d'harmonique ! La valeur de l'angle de rotation a été ramenée dans l'intervalle \( [-\pi, +\pi] \). Les pics que vous observez sur l'évolution du pendule correspondent à des rotations complètes du pendule sur son axe. En bref, il fait n'importe quoi !

C'est confirmé par l'aspect du portrait de phase. Les trajectoires se coupents sur le plan de phase, ce qui est caractéristique du régime chaotique. A vrai dire, le portrait de phase n'est guère exploitable pour analyser la dynamique du système et nous verrons mieux les choses sur sa section de Poincaré.

Spectre de puissance et section de Poincaré

Toujours pour A = 1.2, voici les tracés du spectre de puissance et de la section de Poincaré :

Spectre et SP AF = 1.5

Le spectre de puissance du système est typique des systèmes non linéaires chaotiques. On y observe une raie principale qui correspond à la fréquence de forçage, ce qui le différencie d'un spectre de signal aléatoire et quelques harmoniques de cette raie. Et une multitude de raies d'amplitude variable. On observe que l'amplitude de ces raies est approximativement proportionnelle à l'inverse de leur fréquence. C'est le bruit en "1/f", que l'on déjà rencontré sur TangenteX, à propos du modèle du tas de sable. Ce bruit est typique des systèmes stochastiques. On observe donc ici directement l'apparence aléatoire que peuvent prendre des systèmes déterministes chaotiques.

La section de Poincaré est aussi caractéristique. Comme vous le constatez, les points d'intersection ne sont pas répartis au hasard mais dessinent une forme caractéristique du régime dynamique du système que l'on nomme "attracteur étrange". L'étude des attracteurs étranges d'un système dynamique donne des renseignements très riches sur son comportement. Vous en trouverez des exemples dans les ouvrages [5] à [8] de la bibliographie (voir en fin de page). Remarquons que l'observation de cet attracteur étrange n'est possible que sur la section de Poincaré : sans cet outil, son existence ne serait pas visible.

Exposant de Lyapunov

Comme je l'ai indiqué plus haut, la valeur de l'exposant de Lyapunov est un critère de détermination de la dynamique chaotique d'un système. Pour estimer cette valeur, je vais tracer l'évolution de deux trajectoires dont les conditions initiales seront très proches. Dans l'exemple ci-dessous, la différence entre les deux angles initiaux est de 10-3 radian. Tous les autres paramètres sont égaux, par ailleurs. L'amplitude de forçage est A = 1.5, dans le domaine du régime chaotique comme l'indique le diagramme de bifurcation ci-dessus.

J'ai tracé les deux courbes d'évolution de couleur différente pour observer facilement la divergence entre les deux :

Exposant de Lyapunov - Evolution

Les deux courbes divergent rapidement, au bout de 12 secondes environ (la période de forçage est d'environ 7,7 s), pour devenir complétement différentes après 35 secondes. Si vous choisissiez une valeur plus petite de l'écart des conditions initiales, la divergence prendrait plus de temps, mais apparaîtrait inévitablement. Notre pendule physique est un système très sensible aux conditions initiales.

Il est intéressant de mesurer l'écart entre les deux courbes en fonction du temps. Voici la courbe de ces écarts, figurés sur une échelle logarithmique sur l'axe des ordonnées :

Exposant de Lyapunov - Distance

Les pointes vers le haut que vous observez sont dues au modulo \( 2 \pi \) sur l'angle déjà signalé. Les pointes vers le bas correspondent aux instants où les deux courbes se croisent pour une même valeur de l'angle d'oscillation.

La courbe d'écart croît linéairement entre 0 et 35 secondes environ puis plafonne sur une asymptote, qui correspond à la taille maximum que les courbes peuvent explorer (\( 2 \pi \)). L'échelle des ordonnées étant logarithmique, on en déduit que l'écart entre les deux courbes croît exponentiellement entre 0 et 35 secondes. La mesure de la pente de cette courbe d'écart permet d'estimer expérimentalement la valeur de l'exposant de Lyapunov. On peut aussi l'estimer par le calcul, avec le code suivant :

# calcul du coefficent de Lyapunov
N = distance.size
lya = 0.0
for i in range(N):
    lya += log(abs(distance[i]/epsilon))
lya = lya/N

comme dans le script Lyapunov.py.

Nous constatons qu'il est positif, ce qui confirme, s'il en était besoin, le caractère chaotique de la dynamique de notre pendule physique.

Bibliographie

La dynamique des systèmes non linéaires et le chaos font l'objet d'une très importante littérature. Pour une introduction :

  • [1] "Hasard et Chaos" de David Ruelle.
  • [2] "Chaos et stabilité" de Jérôme Buzzi.
  • [3] "Chaos" de Hans Jürgen Korsch.

et les ouvrages plus techniques :

  • [4] "Des équations différentielles aux systèmes dynamiques" tomes 1 et 2, de Robert Roussarie et Jean Roux.
  • [5] "Differential Equations, Dynamical Systems, and Introduction to Chaos" de Morris W.Hirsch.
  • [6] "Instabilités, Chaos et Turbulence" de Paul Manneville.
  • [7] "Non linear dynamics and chaos" de Steven H.Strogatz.
  • [8] "The nature of chaos" de Tom Mullin.
  • [9] "Dynamical Systems with Applications using Python".

parmi les plus connus...

J'ajouterai à cette liste un ouvrage qui explore le paradigme du pendule en mécanique : "The Pendulum" de Gregory L.Baker. Un vrai bijou, qui dépasse les études théoriques pour aborder les détails pratiques de l'expérimentation avec force schémas. Voir en particulier le schéma et la photo d'un pendule physique dit de "Blackburn-Smith" page 124.

Contenu et design par Dominique Lefebvre - www.tangenteX.com juillet 2020    Licence Creative Commons    Contact :