Lotta per la sopravvivenza

Il ciclo della vita.

La dinamica degli ecosistemi è uno degli argomenti più interessanti, in cui la fisica può essere applicata con successo per studiare l’evoluzione delle popolazioni di prede e predatori.

 

Supponiamo che nel nostro ecosistema ci siano solo due specie in lotta per la sopravvivenza: i pesci (le prede) e gli squali (i predatori).

I pesci hanno a disposizione risorse di cibo praticamente illimitate, se lasciate proliferare si moltiplicano molto rapidamente.

Viceversa gli scquali si cibano solo di pesci, la loro fonte di cibo è quindi molto più scarseggiante, e una sovrappopolazione di predatori causera la loro estinzione per mancanza di cibo. Con una simpatica simulazione è possibile visualizzare la dinamica.

In questa simulazione i quadratini verdi rappresentano i pesci, che si moltiplicano rapidamente, in rosso gli squali, se mangiano pesci si riproducono, altrimenti muoino. Ecco qui un video simpatico che mostra questo ecosistema in azione:

Come si può osservare in numero di pesci e squali aumenta e dinimuisce periodicamente. Se riportiamo su un grafico il numero di prede e quello di predatori notiamo molto bene questo comportamento:

 

Prede e predatori

Prede e predatori

La fisica può descrivere questo tipo di ecosistemi? Si!
Il numero di pesci all’interno del reticolo lo indichiamo con n, il numero di scquali con u.

Scriviamo un’equazione che descriva il comportamento medio di squali e pesci: la probabilità di un pesce di generare un figlio è costante (\alpha), mentre la probabilità che venga mangiata è proporzionale alla probabilità che nei paraggi ci sia uno squalo (\beta).

 
  \displaystyle  \dot n = n\left(\alpha - \beta u\right)
 
Allo stesso modo la probabilità di uno squalo di generare un figlio è tanto maggiore quanto più alta è la probabilità di incrociare un pesce (\gamma), ma hanno una probabilità costante di morire (\delta):

  \displaystyle  \dot u = u\left(\gamma n - \delta\right)
 
Queste equazioni (dette di Lotka-Volterra in onore dei due matematici che le formularono per primi) non hanno soluzione analitica, tuttavia con un cambio di variabili possiamo riscriverle in modo da poterne dare un’interpretazione fisica:
 
  \displaystyle p = \ln n \qquad q = \ln u

  \displaystyle \dot p = \frac{\dot n}{n} \qquad \dot q = \frac{\dot u}{u}
 

Con questa sostituzione il sistema di equazioni che descrive il comportamento medio dei nostri “pesciolini” è:
 
  \displaystyle \left\{\begin{array}{l}  \dot p = \alpha - \beta e^q \\  \dot q = \gamma e^p - \delta  \end{array}\right.
 

Queste due equazioni differenziali adesso descrivono un sistema hamiltoniano canonico:
 
  \displaystyle \left\{\begin{array}{l}  \displaystyle  \dot p = -\frac{\partial H}{\partial q}\\  \\  \displaystyle  \dot q = \frac{\partial H}{\partial p}  \end{array}\right.
 

Si può facilmente ricavare l’hamiltoniana di questo sistema:
 
  \displaystyle  H = -\alpha q -\delta p + \beta e^q + \gamma e^p
 

Quindi è possibile definire un “energia” dell’ecosistema, che è in media conservata (ricordiamoci che le equazioni che abbiamo scritto valgono in media, non tengono conto di piccole fluttuazioni locali della popolazione).
Risostituendo le variabili originali otteniamo che l’energia del sistema é:
 
  \displaystyle  H = \beta u + \gamma n - \alpha \ln u - \delta \ln n
 
Se il numero di pesci e squali è molto alto questo implica che la loro somma, pesata sui coefficienti \beta e \gamma si conserva.

\beta era la probabilità che, data un interazione tra preda e predatore, la preda muoia, mentre \gamma è la probabilità che, data un’interazione tra predatore abbia il sopravvento e quanto questa interazione sia favorevole alla sua procreazione.

Se dividiamo tutto per \beta otteniamo:
  \displaystyle  H' \approx u + \frac{\gamma}{\beta} n

Dove adesso \frac{\gamma}{\beta} è un coefficiente che tiene conto di quanto influenza il cibarsi dei predatori sulla loro crescita. La condizione \gamma = \beta vuol dire che ogni volta che un predatore acchiappa la preda, fa anche un figlio.

Questa condizione tuttavia è molto difficile da raggiungere in pratica, ecco la stima dei parametri per
una delle simulazioni effettuate:
 

  \displaystyle  \omega = \frac{\gamma}{\beta} \qquad \xi =  \frac{\alpha}{\beta} \qquad \zeta = \frac{\delta}{\beta}
 

I valori ottenuti numericamente dalla simulazione sono:
  \displaystyle  \omega \approx 0.29 \qquad \xi \approx 0.36 \qquad \zeta \approx 0.069

Come si può vedere il modello di Lotka-Volterra descrive molto bene l’andamento che è stato misurato.
Conservazione

Oscillatore superarmonico

DomandaCosa hanno in comune un bambino che gioca sull’altalena, una molecola di ossigeno e un matto che fa Bungee Jumping?

Tutti e tre oscillano. La capacità di un sistema fisico di oscillare è una caratteristica molto comune. Ciascuno dei tre esempi risente di forze molto diverse tra loro.
L’altalena del bambino si muove grazie alla forza peso e alla reazione vincolare della catena a cui è sospesa, la molecola di ossigeno vibra a causa del legame coovalente che unisce i due atomi, e il matto che si butta da un ponte con il suo elastico da Bungee Jumping oscilla sotto l’azione combinata della forza di gravità e l’elastico a cui è appeso.

Quando le oscillazioni sono abbastanza piccole, una sola legge fisica accomuna tutti questi moti così diversi tra loro: l’oscillatore armonico. Qualche scienziato afferma addirittura che tutta la fisica si possa ricondurre all’oscillatore armonico.

 

Ma è davvero così?

 

Esaminiamo una forza di richiamo super elastica, di tipo:
   \displaystyle F = -kx^3

Dalla seconda legge della dinamica possiamo ricavare un’equazione che descrive correttamente il moto di un punto materiale soggetto a questa forza:
   \displaystyle F = ma = -kx^3

   \displaystyle m\ddot x = -kx^3

Dove con \ddot x abbiamo indicato la derivata seconda della posizione fatta rispetto al tempo due volte.

Questa equazione differenziale non ammette soluzione analitica. Questa forza è conservativa, e l’energia potenziale è:
   \displaystyle V(x) = \frac{1}{4} k x^4

Il profilo è mostrato in questa figura:

PotenzialeSuperarmonico

 

Se proviamo ad approssimare questo profilo con un potenziale armonico otteniamo una spiacevole sorpresa. Infatti il potenziale superarmonico non è approssimabile con nessun polininomio apparte se stesso. Questo è un esempio di sistema fisico oscillante che non è riconducibile ad un oscillatore armonico!

Non solo non siamo in grado di risolvere l’equazione differenziale che lo definisce, ma per questo potenziale è impossibile tentare un qualsiasi approccio di tipo perturbativo, poiché il più piccolo termine perturbativo non nullo è l’intera forza.

Ci dobbiamo arrendere di fronte a tali difficoltà? Siccome siamo fisici e non matematici, possiamo sfruttare le proprietà fisiche di questo sistema per capirci qualcosa.

Sappiamo che la soluzione di questo potenziale sarà un oscillatore, e vogliamo chiederci quale sia il suo periodo. In generale questo sarà funzione degli unici parametri che compaiono nell’equazione differenziale, e delle condizioni inziali, che qui riassumiamo nella variabile x_0 che rappresenta l’ampiezza di oscillazione:
   \displaystyle T = f(m, k, x_0)

Questa funzione la possiamo spezzare in due parti, una che determinerà la dimensione fisica del periodo (i secondi) e una adimensionale che dipenderà dal dettaglio matematico della soluzione:
   \displaystyle T = m^\alpha k^\beta x_0^\gamma C(m, k, x_0)

Dall’equazione differenziale di partenza possiamo ottenere le dimensioni dei parametri:
   \displaystyle m = [Kg] \qquad k = [\frac{N}{m^3}] = [\frac{Kg}{m^2 s}] \qquad x_0 = [m]

Imponiamo che il periodo di oscillazione si misuri in secondi:

   \displaystyle [s] = [Kg]^\alpha [\frac{Kg}{m^2 s^2}]^\beta [m]^\gamma

Da cui otteniamo il sistema per i coefficienti:

   \displaystyle\left\{\begin{array}{l}   \alpha + \beta = 0 \\   -2\beta + \gamma = 0 \\   -2\beta = 1   \end{array}\right.

Che può essere risolto facilmente:

   \displaystyle \beta = -\frac{1}{2} \qquad \alpha = \frac{1}{2} \qquad \gamma = -1

Da cui abbiamo trovato la dipendenza dimensionale per il periodo:

   \displaystyle T = \frac{1}{x_0}\sqrt{\frac{m}{k}} C(m, k, x_0)

Perché il termine adimensionale sia realmente funzione di quei parametri occorre che ci sia un modo per combinarli in modo da ottenere un risultato adimensionale. Ripetiamo quindi questa prova:

   \displaystyle0 = [Kg]^\alpha [\frac{Kg}{m^2 s^2}]^\beta [m]^\gamma

Da cui otteniamo il sistema:

   \displaystyle\left\{\begin{array}{l}   \alpha + \beta = 0 \\   -2\beta + \gamma = 0 \\   -2\beta = 0   \end{array}\right.

Questo sistema ammette solo la soluzione:
   \displaystyle\alpha = 0 \qquad \beta = 0 \qquad \gamma = 0

Abbiamo dimostrato che il termine C non può dipendere da nessun parametro, è quindi una costante. Da sole considerazioni di tipo dimensionale abbiamo ricavato la formula del periodo dell’oscillatore superarmonico, senza sapere nulla sulla soluzione matematica!

   \displaystyle T = C\frac{1}{x_0}\sqrt{\frac{m}{k}}

Adesso il valore numerico di può essere valutato risolvendo numericamente l’equazione differenziale con qualunque set di parametri iniziale e valutando il periodo numericamente.

Riportiamo in figura la soluzione numerica:

IntegrazioneNumericaDa cui si ottiene una stima di C pari a:

   \displaystyle C \approx 7.42

Dalla figura dell’integrale numerico si vede che rispetto al moto armonico i potenziali superarmonici tendono ad avere una traiettoria più spigolosa, che assomiglia ad un segnale triangolare. Questo perché la potenza superarmonica tende ad appiattire il potenziale (n pari):
   \displaystyle V = \frac{1}{n}kx^{n}\qquad \stackrel{n \rightarrow \infty}{\longrightarrow}\qquad \left\{\begin{array}{lr} \infty & x < -x_0 \mbox{ o } x > x_0 \\ 0 & x \in [-x_0,x_0] \end{array}\right.

Il potenziale assomiglia sempre più ad una buca infinita per n che cresce. All’interno della buca l’oggetto non è più soggetto a forze e tenderà a muoversi di moto rettilineo uniforme, e rimbalzando sulle pareti.

Possiamo calcolare facilmente il periodo per n \rightarrow \infty. Calcoliamo la velocità a cui si muove dentro la buca:
   \displaystyle\frac{1}{2}mv^2 = \frac{1}{n}kx_0^{n}

   \displaystyle v = \sqrt{\frac{2k x_0^n}{n m}}   \qquad   T= \frac{4x_0}{v}

Calcoliamo il periodo usando la legge del moto rettilineo uniforme:

   \displaystyle   T \stackrel{n\rightarrow\infty}{\longrightarrow} 4 \sqrt{\frac{nm}{2k}} x_0^{1 - \frac{n}{2}}

Per il caso che abbiamo trattato (n=4) la costante C approssimata ha un valore di 5.66, neanche troppo distante dal valore reale pari a 7.4. Questo è uno sviluppo corretto per grandi valori di n.

Ad esempio per n = 50 il valore di C stimato da questa relazione è 20.0, quello reale è 20.5. Riportiamo in figura il risultato dell’integrazione per n = 50:

PotenzialeTriangolare

 

Esiste un altro modo interessante per ottenere un espressione analitica di C nel caso n = 4.
Sfruttiamo la conservazione dell’energia meccanica:

   \displaystyle\frac{1}{2}mv^2 + \frac{1}{4}kx^4 = E      \displaystyle v = \sqrt{\frac{2}{m}\left(E - \frac{1}{4}kx^4\right)}

Questa è un altra equazione differenziale. Se ricordiamo che mezzo periodo è il tempo che il nostro oscillatore impiega per spostarsi dalla posizione -x_0 alla posizione x_0 possiamo integrare tutto quanto e ottenere:

   \displaystyle\int_{x_0}^{-x_0} \frac{dx}{\sqrt{\frac{2}{m}\left(E - \frac{1}{4}kx^4\right)}} = \int_0^{\frac{T}{2}}dt      \displaystyle \frac{T}{2} = \sqrt{\frac{m}{2}}\int_{-x_0}^{x_0} \frac{dx}{\sqrt{E - \frac{1}{4}kx^4}}

Ora l’energia totale E del sistema è pari all’energia potenziale superarmonica che l’oscillatore ha nel punto di massima ampiezza:

   \displaystyle E = \frac{1}{4}kx_0^4

 

   \displaystyle T = 2\sqrt{\frac{m}{2}}\frac{2}{x_0^2\sqrt{k}}\int_{-x_0}^{x_0} \frac{dx}{\sqrt{1 - \left(\frac{x}{x_0}\right)^4}}

Per rendere adimensionale l’integrale e ricavare il parametro dobbiamo effettuare il seguente cambiamento di variabili per l’integrale

   \displaystyle \xi = \frac{x}{x_0} \qquad dx = x_0d\xi

 

   \displaystyle T = \frac{1}{x_0}\sqrt{\frac{m}{k}}\sqrt{8} \int_{-1}^1 \frac{d\xi}{\sqrt{1 - \xi^4}}

 

Abbiamo trovato un’espressione analitica per il termine C costante adimensionale nell’espressione del periodo:

   \displaystyle T = \frac{1}{x_0}\sqrt{\frac{m}{k}} C

Questo termine può essere integrato numericamente:

   \displaystyle C = \sqrt{8} \int_{-1}^1 \frac{d\xi}{\sqrt{1 - \xi^4}} \approx 7.42

In accordo con quanto ricavato risolvendo numericamente l’equazione differenziale.