Pendolo semplice, di sicuro ti sarà capitato di vederne uno, giocarci e magari guardare qualche video in cui se ne analizza il comportamento. In questo articolo mi pongo quindi l’obiettivo di provare a formalizzare questa dinamica, analizzando l’equazione differenziale che ne regola il comportamento, disegnando il ritratto di fase e studiandone gli equilibri.
Cercherò di aiutarmi il più possibile con immagini e video perché so che possono aiutare a chiarire le idee, però allo stesso tempo proverò ad usare una procedura rigorosa nell’analisi della dinamica e del ritratto di fase di questo sistema conservativo. Probabilmente in futuro ci sarà tempo per studiare il pendolo smorzato, il pendolo doppio o altre variabili di questo sistema dinamico, ma ora direi che è fondamentale partire dal sistema più semplice, che in realtà richiede comunque un bello studio, ricco di possibili approfondimenti.
L’equazione del pendolo semplice può essere ricavata a partire dal $m\ddot{x}=k\sin{x}$
dove $x$ è l’angolo formato dal pendolo con l’asse verticale e con $\ddot{x}$ intendiamo la derivata seconda di $x$ rispetto al tempo, ovvero l’accelerazione angolare. Può essere che ti sia capitato di vedere l’equazione del pendolo sotto altra forma, ma ho deciso di scriverla con una generica costante $k$ e non coinvolgendo la lunghezza del pendolo, l’accelerazione di gravità e quant’altro per astrarre il più possibile e concentrarci sull’analisi matematica del sistema. Detto ciò, anche nelle forme più “fisiche” di questa equazione, non cambia nulla rispetto a quello che vedremo.
Ah…dimenticavo! Per risolvere un’equazione differenziale chiaramente ci serve una condizione iniziale e, essendo l’equazione di nostro interesse del secondo ordine, ci serve una velocità angolare iniziale ed un angolo iniziale, che chiameremo $(x_0,\omega_0)$. Per tracciare il ritratto di fase lasceremo questi valori liberi, mentre per vedere qualche dinamica particolare, che simulerò con un semplice codice Matlab (che troverai più in basso nell’articolo), ci concentreremo su qualche valore specifico.
Prima di passare all’analisi matematica, penso possa esserti utile una breve introduzione al problema davvero ben fatta da questo video:
Grafico qualitativo del potenziale
Partiamo dal definire cosa intendiamo per sistema conservativo, senza però addentrarci troppo nei dettagli e limitandoci quindi a quello che ci servirà per proseguire nello studio.
Diciamo il sistema definito dalla precedente equazione differenziale essere CONSERVATIVO perché esiste una funzione, $V(x) = k\cos{x}$ nel nostro caso, tale che $\ddot{x}=-V'(x)$, dove $V(x)$ è detto potenziale del sistema.
La nostra equazione differenziale del secondo ordine, è quindi in grado di definire un cosiddetto campo vettoriale del secondo ordine come segue:
$X=X(x,v) = \begin{bmatrix} mv \\ -V'(x) \end{bmatrix}$
A questo punto, per poter tracciare il ritratto di fase, sfruttiamo ulteriormente il fatto che il nostro sistema sia conservativo introducendo un importante integrale primo del nostro campo vettoriale (se non sai cosa sia un integrale primo leggi il mio precedente articolo qui: https://www.mathone.it/integrale-primo/) che è l’energia meccanica, l’energia totale ovvero energia cinetica sommata ad energia potenziale.
Come infatti avrai studiato probabilmente al corso di Fisica 1, ma che comunque mostreremo brevemente qui di seguito, in un sistema conservativo l’energia totale si conserva.
Facciamo un richiamino semplice e definiamo l’energia totale del sistema, verificando poi che esso si conservi effettivamente lungo le curve integrali del campo $X$, ovvero che sia davvero un integrale primo del sistema:
$E=E(x,v) = K(x.v) + V(x) = 1/2mv^2 + V(x)$, di conseguenza possiamo calcolare la derivata di Lie di $E$ lungo il campo $X$ (se non sai cosa sia ti rimando all’articolo citato in precedenza) come segue:
$\mathcal{L}_X(E) = \nabla E \cdot X = \begin{bmatrix} V'(x) \\ mv \end{bmatrix} \cdot \begin{bmatrix} mv \\ -V'(x) \end{bmatrix} = 0 $
Ragioniamo quindi un attimo…cosa sono le orbite di un sistema? Detta semplicemente, per quanto ci riguarda sono curve del piano $x-v$ invarianti rispetto al flusso del campo vettoriale $X$. Sono quindi insiemi 1-dimensionali. Ma siccome la dinamica del sistema definito da $X$ evolve in uno spazio 2-dimensionale (piano $x-v$) e gli insiemi di livello di $E$ sono invarianti, ho una riduzione di un grado di libertà del sistema e gli insiemi di livello di $E$ sono esattamente le orbite del sistema.
Se fossimo quindi in grado, al variare di $c\in\mathbb{R}$, di disegnare gli insiemi
$\Sigma_c = \{ (x,v) : E(x,v) = c\}$ saremmo a posto, perché avremmo disegnato correttamente il ritratto di fase del sistema. Ma perché complicarci la vita?
Infatti è possibile risalire alle stesse orbite solamente disegnando il grafico qualitativo del potenziale $V$, per capire il motivo pensiamoci un po’.
La coppia $(\bar{x},\bar{v})$ appartiene all’insieme $\Sigma_c$ se e solo se $E(\bar{x},\bar{v})=c$. Benissimo, ma è anche sempre vero che $K(x,v) = 1/2 mv^2\geq 0$ per ogni $v$, quindi necessariamente la coppia $(\bar{x},\bar{v})$ per poter stare in $\Sigma_c$ deve essere tale che $\bar{x}\in \sigma_c = \{ x : V(x)\leq c\}$, altrimenti è impossibile ottenere $c$ come somma del potenziale con l’energia cinetica. Ecco quindi che ottieniamo quanto anticipato prima, infatti da ciò segue che per disegnare $\Sigma_c$ basta capire com’è fatto $\sigma_c$ che è molto più semplice da studiare visto che è il sottolivello $c$ di una funzione di una sola variabile.
Ecco quindi il grafico qualitativo del potenziale del nostro sistema qui di seguito:
Per tracciare quindi le varie orbite, corrispondenti ai diversi livelli energetici, nel prossimo paragrafo analizzeremo le intersezioni delle rette di forma $y=c$ con questo grafico.
Da grafico qualitativo a ritratto di fase
Come abbiamo già introdotto in precedenza, per ottenere il ritratto di fase del sistema in analisi, è sufficiente studiare i sottoinsiemi di livello della funzione potenziale. Ho quindi evidenziato qui di seguito i vari livelli energetici interessanti, rappresentandoli con 5 colori diversi.
Per proseguire, dobbiamo utilizzare alcune regole utili per passare dal grafico qualitativo del potenziale ad una specifica orbita o insieme di orbite. Volendo si potrebbero ricavare anche con un’interpretazione fisica queste “regole” ma per il momento mi limito ad elencare quelle che ci interessano, per poi applicarle al nostro esempio.
Intanto definiamo sottoinsieme di livello, in corrispondenza del livello energetico $E=c$ il seguente insieme $\sigma_c = \{x\in\mathbb{R}: V(x)\leq c\}$. In base alle proprietà di questo insieme, prevalentemente proprietà topologiche, si distinguono le varie tipologie di orbite. Inoltre, nel caso tale insieme non sia connesso (se non sai cosa significhi essere connesso vai a leggerti una spiegazione qui), si può procedere all’analisi sulle sue componesse, ovvero sui suoi sottoinsiemi connessi e disgiunti, la cui unione restituisce l’insieme disconnesso di partenza.
Ecco quindi le “regole” che andremo a seguire:
- Se $\sigma_c = [a,b]$ è limitato, allora ad esso è associata un orbita periodica, simmetrica rispetto all’asse x, che la interseca con una tangente verticale. Vediamo un esempio:
Dove i punti di intersezione con l’asse $x$ corrispondono alla $a$ e $b$ del grafico di sinistra.
- Se $\sigma_c=[-\infty,+\infty]$ allora l’orbita associata è illimitata, per tempi infinitamente grandi e piccoli questa andrà rispettivamente come $x\rightarrow \pm\infty$ ed è simmetrica rispetto all’asse $x$, senza intersecarla, come puoi vedere di seguito:
- Nel caso della retta verde chiaro, abbiamo ancora che $\sigma_c=\mathbb{R}$ ma, inoltre, abbiamo un unico punto in cui $V(x) = c$, che è un equilibrio instabile. In questo caso la situazione è un po’ più complicata e avremo un punto chiamato di equilibrio instabile, è un punto iperbolico e in sua corrispondenza abbiamo 5 orbite, ottenendo quanto segue:
- Nel caso in cui il sottoinsieme di livello abbia un solo punto, come con la retta viola che puoi vedere di seguito, avremo un equilibrio stabile e l’unica orbita ad esso associata è il punto stesso.
- Per concludere, nel caso in cui il sottoinsieme di livello sia vuoto, come presentato qui sotto, avremo che non c’è nessun’orbita associata.
Benissimo, ora abbiamo tutti gli strumenti per “unire i puntini” e passare dal grafico della funzione potenziale al ritratto di fase. Per farlo basta semplicemente assemblare i 5 pezzi visti qui sopra e accorgersi della periodicità del nostro potenziale. Infatti è sufficiente analizzare il ritratto di fase tra due equilibri instabili e poi ripetere per periodicità.
Vediamo quindi cosa si ottiene assemblando il tutto! Il risultato che trovi qui sotto l’ho ottenuto con un servizio gratuito online che permette di rappresentare ritratti di fase di campi vettoriali, lo trovi qui: https://aeb019.hosted.uark.edu/pplane.html . Più in basso troverai invece una video simulazione del codice che ho scritto in Matlab per ottenere lo stesso risultato.
Risoluzione esatta equazione del pendolo
Una volta letto il precedente paragrafo, sono certo che ti sarai convinto che possiamo tranquillamente scegliere una condizione iniziale a caso e sapremo quale sarà l’orbita ad essa associata e quindi la dinamica del sistema.
Un altro modo per ottenere lo stesso risultato è risolvere direttamente il problema di Cauchy seguente:
$m\ddot{x}=k\sin{x}$ con $x(0)=x_0$ e $\dot{x}(0)=\omega_0$
L’equazione differenziale del secondo ordine in analisi la si può risolvere in vari modi, probabilmente il più classico è il metodo di somiglianza, anche noto come metodo degli annichilatori, che trovi spiegato nei video qui di seguito (nel primo è introdotto il metodo, mentre nel secondo si tratta l’esempio che ci interessa):
Altro modo per risolverla potrebbe essere ricondursi ad un sistema del primo ordine e lavorare di conseguenza. Per evitare di dilungarmi eccessivamente, se ti interessa scrivi pure un commento a questo articolo che ti mostrerò anche questa procedura.
Otteniamo quindi $x(t) = -k/m (\sin{t}-t(v_0 m/k+1))+x_0$ e quindi il suo grafico potrà descrivere tranquillamente la dinamica del sistema, in alternativa a quanto fornito dal ritratto di fase.
Qui ho semplicemente riportato il grafico della soluzione per una particolare scelta della condizione iniziale $(x_0,v_0)=(0,1)$, come puoi vedere nel grafico. Infatti questa curva passa per l’origine e la tangente in essa è a pendenza 1, proprio come $v_0$. Ci tenevo a riportare questo grafico per chiarire che nonostante il comportamento del sistema dinamico che ci interessa sia periodico, la fase o comunque la variabile $x=x(t)$ continua a crescere nel tempo.
Implementazione numerica dell’equazione del pendolo
<?matlab
clear all
close all
k = 1;
m = 1;
f = @(t,y) [y(2);
k/m * sin(y(1))];
t0 = 0; T = 8;
n = 10001;
timespan = linspace(t0,T,n);
dt = timespan(2)-timespan(1);
initial = [pi .2 pi/2 pi-.1 1 0.1;
0 1 .5 1 3 0.1 ];
for j = 1:length(initial(1,:))
x0 = initial(1,j);
v0 = initial(2,j);
y0 = [x0;v0];
plot(x0,v0,'ko','MarkerSize',5)
t = t0;
y = zeros(2,length(timespan));
y(:,1) = y0;
for i = 1:length(timespan)-1
y(:,i+1) = y(:,i) + dt * f(t,y(:,i));
t = t + dt;
end
if j>1
plot(y(1,:),y(2,:),'r*','MarkerSize',1);
else
plot(y(1,:),y(2,:),'ro','MarkerSize',6);
end
hold on;
end
?>
[/matlab]E qui ho riportato l’esito di questo codice, in cui vengono simulate 6 orbite, partendo da differenti condizioni iniziali. Chiaramente non essendo curato e ottimizzato gli output sono molto meno chiari rispetto a quelli riportati nelle sezioni precedenti. Tuttavia ho preferito aggiungere questo paragrafo per mostrare, anche a chi non ha mai scritto nessuna implementazione numerica, come sia possibile “simulare” l’andamento di un particolare sistema dinamico.
Giusto perché si sappia, per integrare nel tempo ho usato Eulero Esplicito, non è senz’altro la miglior strada percorribile ma è la più semplice, da cui iniziare per ottenere qualche risultato tangibile.
Lascia un commento