Poincaré Maps¶
da PowerShell -- jupyter nbconvert --to html notebook.ipynb
Mappe di Poincaré e doppio pendolo¶
Motivazione
Lo studio qualitativo dei sistemi dinamici non lineari è spesso ostacolato dall’elevata dimensionalità dello spazio delle fasi. Nel caso dei sistemi hamiltoniani a due gradi di libertà , come il doppio pendolo piano, lo spazio delle fasi è quadridimensionale, rendendo impossibile una visualizzazione diretta delle traiettorie. La $\textit{mappa di Poincaré}$ fornisce uno strumento fondamentale per ridurre la dimensionalità del problema e analizzarne la struttura dinamica.
Il doppio pendolo come sistema hamiltoniano
Il doppio pendolo piano è costituito da due masse puntiformi $ m_1, m_2 $ collegate da aste rigide e prive di massa di lunghezza $ l_1, l_2 $, vincolate a muoversi in un piano verticale sotto l’azione della gravità . Il sistema è descritto da due coordinate generalizzate, gli angoli $ \theta_1, \theta_2 $, e dalle rispettive velocità angolari $ \dot{\theta}_1, \dot{\theta}_2 $.
Lo stato del sistema è quindi rappresentato dal vettore $$\large (\theta_1, \dot{\theta}_1, \theta_2, \dot{\theta}_2) \in \mathbb{R}^4, $$ che evolve secondo equazioni del moto non lineari derivate dal formalismo lagrangiano. Il sistema è conservativo e possiede un integrale del moto, l’energia totale $ E $, ma non è integrabile in senso generale.
Definizione della mappa di Poincaré
Sia dato un sistema dinamico continuo $$\large \dot{\mathbf{x}} = \mathbf{f}(\mathbf{x}), \qquad \mathbf{x} \in \mathbb{R}^n. $$ Una mappa di Poincaré è definita scegliendo una superficie di sezione $ \Sigma \subset \mathbb{R}^n $ di codimensione uno, trasversale al flusso. La mappa associa a un punto $ \mathbf{x}_k \in \Sigma $ il punto $ \mathbf{x}_{k+1} \in \Sigma $ ottenuto alla successiva intersezione della traiettoria con $\Sigma $.
Nel caso del doppio pendolo, una scelta naturale della sezione è $$\large \Sigma = \{ \theta_2 = 0,\ \dot{\theta}_2 > 0 \}, $$ cioè l’insieme degli stati in cui il secondo pendolo attraversa la verticale con velocità angolare positiva.
Riduzione dimensionale
L’intersezione della traiettoria con la sezione $ \Sigma $ produce una successione discreta di punti nello spazio delle fasi. Poiché l’energia totale è costante, la dinamica è confinata a una varietà tridimensionale, e la mappa di Poincaré risulta bidimensionale.
In pratica, si rappresentano le coppie $$\large (\theta_1, \dot{\theta}_1) $$ valutate a ogni attraversamento della sezione. La mappa di Poincaré è quindi una trasformazione discreta del piano $(\theta_1, \dot{\theta}_1)$.
Interpretazione geometrica
La struttura della mappa di Poincaré riflette direttamente la natura del moto:
- $\textbf{Curve chiuse lisce}$: moto quasi-periodico, associato a tori invarianti nello spazio delle fasi.
- $\textbf{Catene di isole}$: risonanze tra le frequenze fondamentali del sistema.
- $\textbf{Regioni riempite in modo irregolare}$: moto caotico, caratterizzato da sensibilità alle condizioni iniziali.
- $\textbf{Curve spezzate o parzialmente riempite}$: regime misto, con sopravvivenza parziale dei tori KAM
e presenza di barriere frattali (cantori).
Ruolo dell’energia
Nel doppio pendolo il comportamento dinamico dipende in modo critico dall’energia totale:
- A bassa energia il sistema è quasi integrabile e la mappa di Poincaré mostra curve regolari.
- A energia intermedia compaiono risonanze, isole di stabilità e regioni caotiche locali.
- A energia elevata il caos diventa globale e la sezione viene riempita densamente.
La mappa di Poincaré consente quindi di visualizzare in modo diretto la transizione da moto regolare a moto caotico in un sistema hamiltoniano.
Significato fisico
La mappa di Poincaré non è solo uno strumento grafico, ma rappresenta una riduzione dinamica del problema continuo a uno discreto, rendendo visibili fenomeni fondamentali come: la rottura dei tori invarianti, la diffusione caotica, e la struttura frattale dello spazio delle fasi. Nel caso del doppio pendolo, essa costituisce uno degli esempi paradigmatici di caos deterministico nella meccanica classica.
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import solve_ivp
g = 9.81
l1 = 1.0
l2 = 0.5
m1 = 1.0
m2 = 0.2
def double_pendulum(t, y):
th1, w1, th2, w2 = y
d = th1 - th2
den1 = (m1 + m2)*l1 - m2*l1*np.cos(d)**2
den2 = (l2/l1)*den1
dw1 = (
m2*l1*w1**2*np.sin(d)*np.cos(d)
+ m2*g*np.sin(th2)*np.cos(d)
+ m2*l2*w2**2*np.sin(d)
- (m1 + m2)*g*np.sin(th1)
) / den1
dw2 = (
-m2*l2*w2**2*np.sin(d)*np.cos(d)
+ (m1 + m2)*(g*np.sin(th1)*np.cos(d)
- l1*w1**2*np.sin(d)
- g*np.sin(th2))
) / den2
return [w1, dw1, w2, dw2]
# bassa energia - assenza di caos
y0 = [
0.3, # theta1
0.0, # theta1_dot
0.2, # theta2
0.0 # theta2_dot
]
# solo qulche risonanza è superata - caos confinato con isole
#y0 = [0.5, 0.0, 0.3, 0.0]
# vario di pochissimo, il grafico cambia drasticamente - caos confinato con isole
y0 = [0.499, 0.0, 0.3, 0.0]
# alta energia - solo caos confinato
#y0 = [0.7, 0.0, 0.5, 0.0]
t_max = 400
dt = 0.01
t_eval = np.arange(0, t_max, dt)
sol = solve_ivp(
double_pendulum,
(0, t_max),
y0,
t_eval=t_eval,
rtol=1e-9,
atol=1e-9
)
theta1 = sol.y[0]
z1 = sol.y[1]
theta2 = sol.y[2]
z2 = sol.y[3]
theta1 = (theta1 + np.pi) % (2*np.pi) - np.pi
theta2 = (theta2 + np.pi) % (2*np.pi) - np.pi
theta1_p = []
z1_p = []
for i in range(len(theta2)-1):
if theta2[i] < 0 and theta2[i+1] > 0 and z2[i] > 0:
# interpolazione lineare
a = -theta2[i] / (theta2[i+1] - theta2[i])
th1 = theta1[i] + a*(theta1[i+1] - theta1[i])
w1 = z1[i] + a*(z1[i+1] - z1[i])
theta1_p.append(th1)
z1_p.append(w1)
theta1_p = np.array(theta1_p)
z1_p = np.array(z1_p)
plt.figure(figsize=(6,6))
plt.scatter(theta1_p, z1_p, s=6)
plt.xlabel(r'$\theta_1$')
plt.ylabel(r'$\dot{\theta}_1$')
plt.title(r'Mappa di Poincaré – $\theta_2 = 0,\ \dot{\theta}_2 > 0$')
plt.grid(True)
plt.show()
Rotazione caotica di Hyperion¶
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.animation import FuncAnimation
from IPython.display import HTML
# Costanti
e = 0.1
Iratio = 0.264
a = 1 / np.sqrt(1 - e**2)
b = 1.5 * a * Iratio
pi = np.pi
# Parametri
slices = 20 # numero di sezioni
dt = 2 * pi / slices
orbits = 30
itns = 100 # periodi di Keplero
def r(eta):
"""Orbita kepleriana"""
return (1 - e**2) / (1 + e * np.cos(eta))
def fHyp(z):
"""
RHS equazioni rotazionali
z = [u, v, eta]
"""
u, v, eta = z
return np.array([
a * r(eta)**2 * v,
-b * np.sin(2*u - 2*eta),
1.0
])
# integrazione rk4 - passo singolo
def rk4_step(f, dt, z):
k1 = f(z)
k2 = f(z + 0.5 * dt * k1)
k3 = f(z + 0.5 * dt * k2)
k4 = f(z + dt * k3)
return z + (dt/6) * (k1 + 2*k2 + 2*k3 + k4)
# condizioni iniziali
initz = np.array([
[0, 0, 0],
[0, 1/7, 0],
[0, 9/28, 0],
[0, 3/7, 0],
[0, 15/28, 0],
[0, 6/7, 0],
[0, 8/7, 0],
[0, 10/7, 0],
[0, 12/7, 0],
[0, 2, 0],
[0, 16/7, 0],
[0, 5/2, 0],
[0, 37/14, 0],
[0, 39/14, 0],
[0, 81/28, 0],
[0, 3, 0],
[pi/2, 1/7, 0],
[pi/2, 2/7, 0],
[pi/2, 3/7, 0],
[pi/2, 53/28, 0],
[pi/2, 57/28, 0],
[pi/2, 29/14, 0],
[pi/2, 18/7, 0],
[pi/2, 20/7, 0],
[pi, 9/28, 0],
[pi, 3/7, 0],
[pi, 15/28, 0],
[pi, 6/7, 0],
[pi, 16/7, 0],
[pi, 37/14, 0]
])
selected_slices = [0, slices//3, 2*slices//3]
tbl = {}
maps = [[] for _ in range(slices)]
for n in range(orbits):
z = initz[n].copy()
for i in range(itns):
for j in range(slices):
z = rk4_step(fHyp, dt, z)
maps[j].append([np.mod(z[0], 2*pi), z[1]])
maps = [np.array(m) for m in maps]
fig, ax = plt.subplots(figsize=(5,5))
scat = ax.scatter([], [], s=4, c='black')
ax.set_xlim(0, 2*pi)
ax.set_ylim(0, 3)
ax.set_aspect('equal')
ax.axis('off')
def init():
scat.set_offsets(np.empty((0,2)))
return scat,
def update(frame):
scat.set_offsets(maps[frame])
ax.set_title(f"Sezione {frame+1}/{slices}", fontsize=10)
return scat,
anim = FuncAnimation(
fig,
update,
frames=slices,
init_func=init,
interval=400, # ms tra i frame
blit=True
)
HTML(anim.to_jshtml())
#plt.close(fig)
