Nuclei di memoria in modelli biologici di Volterradigilander.libero.it/matteolesinigo/Matteo/Tesi di...

POLITECNICO DI MILANOFacolta di Ingegneria dei Sistemi

Dipartimento di Matematica F.BrioschiCorso di studi in Ingegneria Matematica

Tesi di Laurea di I livello

Nuclei di memoriain modelli biologici di Volterra

Relatrice: Professoressa Conti Monica

Lesinigo Matteo Zanzottera AnnaMatricola 661380 Matricola 663765

Anno Accademico 2005/2006

Un grande, profondo e sincero grazie alla Professoressa Monica Conti che con la sua presenzacostante, la sua gentilezza e cordialita e le sue correzioni sempre costruttive ci ha aiutato nellastesura di questo elaborato e soprattutto ha arricchito il nostro bagaglio culturale e umano.

Un enorme grazie ai miei genitori che mi hanno sempre sostenuto e incoraggiato durante tutta lamia vita e mi hanno consentito di frequentare questo splendido corso di laurea,Un grazie speciale alla mia fidanzata Vera che da due anni e mezzo mi sostiene con tanta pazienzae continua a farlo,Un grazie al mio piu caro amico Gianpaolo per tutto quello che fa e ha fatto per me in questi anni,Un grazie a tutti i miei parenti che mi vogliono bene e mi incoraggiano con particolare riferimentoalle mie due nonne, Elide e Luciana,Un grazie a tutti i miei amici, in rigoroso ordine alfabetico, Carlo, Claudia, Elisa, Elisabetta,Emiliano, Francesca, Laura e Stefano,E naturalmente un grande grazie alla mia amica Anna con cui ho lavorato benissimo in questi mesida tesisti...

Matteo Lesinigo

Un grazie particolarea Matteo che con la sua determinazione e originalita, mi ha aiutato e continua ad aiutarmi acercare il “modello” adatto ad affrontare i problemi non solo matematicia Mamma, Papa e Paolo che, con i loro consigli, la loro pazienza e il loro affetto mi aiutano a darestabilita a tale “modello”a tutte le amiche e gli amici che tentano invece di destabilizzarlo ma che in realta lo arricchisconocon il loro entusiasmo, la loro esuberanza e le loro critiche.

Anna Zanzottera

Indice

Prefazione iii

1 Modelli matematico-biologici di evoluzione delle popolazioni 11.1 Il modello di Malthus . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1

1.1.1 Introduzione . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.1.2 Analisi della soluzione . . . . . . . . . . . . . . . . . . . . . . . . . . 2

1.2 Il modello Logistico . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21.2.1 Introduzione . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21.2.2 Analisi della soluzione . . . . . . . . . . . . . . . . . . . . . . . . . . 4

1.3 Il modello con memoria per singola popolazione:l’equazione di Verhulst-Volterra . . . . . . . . . . . . . . . . . . . . . . . . . 41.3.1 Introduzione . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41.3.2 Considerazioni sull’equazione di Verhulst-Volterra . . . . . . . . . . 6

1.4 Il sistema di Lotka-Volterra . . . . . . . . . . . . . . . . . . . . . . . . . . . 61.4.1 Introduzione . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61.4.2 Analisi delle soluzioni . . . . . . . . . . . . . . . . . . . . . . . . . . 71.4.3 Teoremi sul sistema di Lotka-Volterra . . . . . . . . . . . . . . . . . 9

1.5 Il sistema di Lotka-Volterra con controllo logistico . . . . . . . . . . . . . . 161.5.1 Introduzione . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 161.5.2 Analisi delle soluzioni . . . . . . . . . . . . . . . . . . . . . . . . . . 16

1.6 Il sistema di Lotka-Volterra con memoria . . . . . . . . . . . . . . . . . . . 171.6.1 Introduzione . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 171.6.2 Considerazioni sul sistema di Lotka-Volterra con memoria . . . . . . 21

2 Nuclei a decadimento esponenziale 232.1 Introduzione . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 232.2 Il nucleo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 242.3 Da sistemi integro-differenziali a sistemi autonomi . . . . . . . . . . . . . . 252.4 Convergenza del nucleo alla delta di Dirac nel senso delle distribuzioni

quando r tende a 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 282.5 Esempi di convoluzioni con funzioni note . . . . . . . . . . . . . . . . . . . . 30

3 L’equazione di Verhulst-Volterra 373.1 Dall’equazione al sistema . . . . . . . . . . . . . . . . . . . . . . . . . . . . 373.2 Analisi degli equilibri del sistema (3.1.3) . . . . . . . . . . . . . . . . . . . . 38

i

ii INDICE

3.3 Confronto con il Teorema 1.3.1 . . . . . . . . . . . . . . . . . . . . . . . . . 423.4 Soluzioni al variare del PDN . . . . . . . . . . . . . . . . . . . . . . . . . . . 49

4 Il sistema di Lotka-Volterra con memoria 554.1 Dal sistema integro-differenziale a quello autonomo . . . . . . . . . . . . . . 554.2 Analisi degli equilibri del sistema (4.1.3) . . . . . . . . . . . . . . . . . . . . 57

4.2.1 Analisi dell’equilibrio nullo (4.2.1a) . . . . . . . . . . . . . . . . . . . 574.2.2 Analisi dell’equilibrio (4.2.1b) (estinzione dei predatori) . . . . . . . 584.2.3 Analisi dell’equilibrio (4.2.1d)(coesistenza di entrambe le specie) . . 604.2.4 Un esempio modellisticamente significativo . . . . . . . . . . . . . . 65

4.3 Limitatezza delle soluzioni . . . . . . . . . . . . . . . . . . . . . . . . . . . . 694.4 Cicli limite . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 734.5 Soluzioni al variare dei PDN . . . . . . . . . . . . . . . . . . . . . . . . . . . 78

4.5.1 Convergenza del sistema a quello senza memoria quando il PDNtende a zero . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78

4.5.2 Convergenza del sistema a quello senza memoria quando il PDNtende a infinito . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79

Conclusioni 85

A Risultati utilizzati nel testo 87A.1 Algebra Lineare . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87A.2 Analisi Funzionale . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 88

A.2.1 Metodi di Linearizzazione per sistemi di EDO non lineari autonomi 89

Elenco delle Figure 93

Bibliografia 95

Prefazione

Il lavoro si propone di esaminare gli aspetti salienti dei modelli matematico-biologici pro-posti da Volterra1 per descrivere l’evoluzione del numero di individui costituenti una po-polazione. Dopo una breve introduzione, il cui scopo e quello di familiarizzare il lettorecon i modelli piu comuni di equazioni differenziali atti a descrivere la dinamica delle popo-lazioni, verranno introdotti i modelli di equazioni integro-differenziali di Verhulst-Volterrae di Lotka-Volterra, cardini del nostro lavoro. Tali modelli verranno studiati nel detta-glio nel caso di nuclei a decadimento esponenziale proponendo numerosi esempi numericiottenuti con il software Mathematica2. Senza troppe pretese la nostra trattazione parteda una accurata analisi degli andamenti qualitativi delle soluzioni numeriche dei modelliintegro-differenziali che, come il lettore avra modo di vedere, risultano di profondo inte-resse modellistico. Qualora possibile, data l’assenza di una teoria ben sviluppata in questocampo (soprattutto per quanto riguarda i sistemi di equazioni integro-differenziali) verran-no proposti risultati teorici o congetture che risultino supportate da evidenza numerica.In una sezione a parte verra spiegata in modo esaustivo la scelta di utilizzare i nucleia decadimento esponenziale. Risultati di riferimento per l’analisi teorica di equazioni osistemi differenziali ordinari sono presenti in appendice.

1Vito Volterra (1860 - 1940) matematico e fisico italiano, uno dei principali fondatori dell’analisi fun-zionale e della connessa teoria delle equazioni integrali, noto anche per i suoi contributi alla biologiamatematica.

2Mathematica 5.2 un marchio registrato della Wolfram Research Inc. Per maggiori dettagli si visiti ilsito www.wolfram.com

iii

Capitolo 1

Modelli matematico-biologici dievoluzione delle popolazioni

Scopo di questo Capitolo e passare in rassegna i modelli matematico-biologici atti a descri-vere l’evoluzione di una popolazione o di due popolazioni interagenti. I modelli vengonopresentati in ordine di complessita. Tale ordine puo anche essere considerato l’ordine cro-nologico con il quale tali modelli sono stati proposti. Per ogni modello vengono sottolineatigli aspetti salienti e si mette in luce il significato biologico dei parametri da cui il modellostesso dipende. Si provvede inoltre a fornire descrizioni di eventuali punti di equilibrio ead integrare la trattazione con alcuni grafici che mettano in luce i comportamenti fonda-mentali dei singoli modelli. Visto che l’interesse della trattazione e strettamente legato alsignificato biologico delle variabili, tutte le considerazioni proposte sottintendono che levariabili che modellizzano il numero di individui siano sempre non negative. I primi mo-delli sono relativi ad una singola popolazione mentre gli ultimi riguardano due popolazioniinteragenti tra loro (prede e predatori). La descrizione dei modelli piu semplici relativi al-l’uno e all’altro caso e seguita dalla deduzione dei modelli con memoria, oggetto principaledella tesi. Questi ultimi verranno ampiamente trattati nel caso di memoria evanescente adecadimento esponenziale nei Capitoli 3 e 4.

1.1 Il modello di Malthus

1.1.1 Introduzione

Nel 1798 Thomas Robert Malthus, economista e storico inglese, propose un primo semplicemodello matematico per descrivere l’evoluzione di una popolazione. Egli considero unasituazione ideale in cui non vi fosse alcun limite alle risorse di cibo e di grandezza dell’ha-bitat per una popolazione di individui. Sotto tali ipotesi Malthus ritenne ragionevole chela variazione del numero di individui in un dato istante temporale t fosse proporzionale alnumero di individui presenti al tempo t medesimo mediante un parametro k che prendeil nome di potenziale biologico. Tale potenziale biologico ha il significato di differenza trail tasso di natalita e quello di mortalita per individuo nell’unita di tempo. Se dunquesi indica con x(t) il numero di individui al tempo t e si suppone che x sia una funzione

1

2 Modelli matematico-biologici di evoluzione delle popolazioni

continua1 e differenziabile si puo giungere a scrivere il modello:

x′(t) = k x(t). (1.1.1)

Il modello di Malthus risulta quindi un’equazione differenziale del primo ordine, lineare acoefficiente costante. Tale equazione, associata al dato iniziale x(t0) = x0 ammette una euna sola soluzione dipendente dal parametro k:

x(t) = x0 exp[k (t− t0)]. (1.1.2)

1.1.2 Analisi della soluzione

Per semplicita considereremo t0 = 0. Questa semplificazione non e assolutamente restrit-tiva visto che l’equazione e tempo invariante2. Si nota subito che esiste un unico equilibriocorrispondente all’origine:

xeq = 0.

Tale equilibrio risulta:

• Globalmente Asintoticamente Stabile se il potenziale biologico risulta negativo. Ilsignificato di tale affermazione e immediato. Nelle ipotesi fatte se il numero diindividui che muore e maggiore del numero di individui che nasce per unita di tempola popolazione tende ad estinguersi.

• Stabile se il potenziale biologico e nullo. Cio significa che il numero di individuisi mantiene limitato (nel caso specifico costante e pari a x0 se il tasso di natalitauguaglia quello di mortalita.

• Instabile se il potenziale biologico e positivo. Cio significa che se nascono piu indi-vidui di quanti ne muoiano la popolazione, avendo a disposizione risorse illimitatecresce a dismisura (in particolare in modo esponenziale).

Gli andamenti fondamentali che caratterizzano il modello evoluzionistico di Malthus sonoriportati per un dato iniziale particolare (x0 = 10) in Figura 1.1.

1.2 Il modello Logistico

1.2.1 Introduzione

Il modello di Malthus risulta irrealistico nella maggior parte dei casi. In effetti l’ipotesidi risorse illimitate sembra poco plausibile e risulta evidente come al crescere del numero

1N.B. Tale approssimazione puo ritenersi ragionevole se il numero di individui e sufficientemente grande.Infatti l’errore relativo che si commette approssimando un valore reale di x(t) all’intero piu vicino (ilnumero di individui nella realta e logicamente un intero) diminuisce al crescere della popolazione. Perpopolazioni piccole sarebbe preferibile optare per equazioni alle ricorrenze o per sistemi dinamici discreti.Gli esempi proposti nel nostro lavoro, pur facendo riferimento sempre ad equazioni differenziali, non semprecontemplano l’analisi di popolazioni grandi. In ogni caso i valori dei parametri e dei dati delle simulazioninumeriche sono sempre stati scelti in modo tale da mettere in luce tutte le caratteristiche principali delmodello in esame.

2Se si cambia l’istante iniziale la soluzione rimane la medesima a meno di traslazioni temporali

Il modello Logistico 3

Figura 1.1: Soluzioni dell’equazione (1.1.1) relative a diversi valori del potenziale biologico k incorrispondenza del dato iniziale x0 = 10. In ordinata viene riportato il numero di individui dellapopolazione mentre in ascissa il tempo. In rosso la soluzione per k = 0.2, in verde per k = 0 e inblu per k = −0.5.

di individui le risorse pro-capite diminuiscano sensibilmente. Tale effetto non e tenutoin considerazione dal modello visto in precedenza e appare logico tentare di modificare ilmodello stesso cosicche incorpori le restrizioni imposte dalla limitatezza delle risorse. Nel1845 Verhulst3 propose un modello nel quale ipotizzo una decrescita lineare in funzionedel numero di individui del tasso relativo di crescita della popolazione. Con questa ipotesila legge di evoluzione (1.1.1) si puo dunque riscrivere nel seguente modo:

x′(t) = k x(t)(1− x(t)

ε

). (1.2.1)

A tale equazione andra naturalmente associata una opportuna condizione iniziale x(0) =x0. L’equazione (1.2.1) e un’equazione differenziale ordinaria del primo ordine non linearea coefficienti costanti. Il parametro k continua ad avere il significato di potenziale biologicodescritto nella sezione 1.1. A tale parametro se ne aggiunge un secondo pari a −k

ε cheha il compito di contrastare la crescita della popolazione. Tale termine che tiene contodella limitatezza delle risorse viene comunemente indicato in letteratura come controllologistico. Nonostante l’equazione sia non lineare ammette una soluzione in forma chiusadata dalla seguente espressione:

x(t) =ε x0 exp[k t]

ε− x0 + x0 exp[k t]. (1.2.2)

Si suppone che sia il potenziale biologico che il parametro ε siano positivi. Cio escludela situazione, del resto assolutamente assurda dal punto di vista modellistico, in cui unapopolazione con potenziale biologico negativo possa crescere a dismisura come previstodalla (1.2.2) nel caso k < 0 e x0 > ε > 0.

3Pierre Francois Verhulst (1804-1849), matematico belga.


Figura 1.2: Soluzioni dell’equazione (1.2.1) corrispondenti a diversi valori del dato iniziale x0. Intutti e tre i casi si e supposto k = 0.5 e ε = 20. In ordinata viene riportato il numero di individuidella popolazione mentre in ascissa il tempo. In rosso la soluzione per x0 = 10, in verde per x0 = 20e in blu per x0 = 30.

1.2.2 Analisi della soluzione

Per prima cosa osserviamo come ora vi siano due punti di equilibrio distinti. Tali punti diequilibrio sono:

xeq1 = 0, (1.2.3a)

xeq2 = ε. (1.2.3b)

Linearizzando4 il sistema e tenendo conto della positivita dei parametri si deduce fa-cilmente che l’equilibrio (1.2.3a) e sempre instabile mentre l’equilibrio (1.2.3b) e sempreglobalmente asintoticamente stabile. A seguito di queste considerazioni il parametro ε puoessere chiamato capacita dell’ambiente. Riassumiamo quanto detto nel seguente risultato:

Microteorema 1.2.1. Qualora il numero di individui iniziale sia strettamente positivo lapopolazione non si estingue e tende alla capacita dell’ambiente.

In Figura 1.2 viene mostrato l’andamento delle soluzioni al variare del dato iniziale incorrispondenza di k = 0.5 e ε = 20.

1.3 Il modello con memoria per singola popolazione:l’equazione di Verhulst-Volterra

1.3.1 Introduzione

I modelli 1.1 e 1.2 possono naturalmente essere raffinati aggiungendo termini di ordinesuperiore al secondo ottenendo varie tipologie di comportamento a seconda dell’ordinemedesimo. Risulta tuttavia di particolare interesse un altro approccio. Considerando

4Per maggiori dettagli circa il metodo di linearizzazione si veda l’appendice A.2.1

l’equazione di Verhulst Volterra 5

i modelli 1.1 e 1.2 si puo anche dire che la derivata logaritmica del numero di individuirisulta nel primo caso esattamente uguale al potenziale biologico k e nel secondo caso ugualealla differenza tra il potenziale biologico stesso e un termine proporzionale al numero diindividui mediante il parametro di controllo logistico. Entrambe queste relazioni sonotempo invarianti e esprimono una correlazione istantanea tra i due membri. L’ipotesidi relazione istantanea tra la derivata logaritmica del numero di individui e i terminievidenziati sembra pero poco realistica. Un esempio banale puo aiutare a comprenderemeglio la necessita di modificare questo tipo di relazione.

Esempio 1.3.1. Consideriamo i tempi di gestazione: la variazione del numero di individuial tempo t e correlata alla quantita di individui al tempo del concepimento e non del parto.Modifichiamo l’equazione (1.1.1) introducendo un ritardo. Le strade praticabili sono varie.Scegliamo di legare il tasso relativo di crescita del numero di individui al tempo t al numerodi individui al tempo t − τ (ove τ > 0 e un opportuno tempo di ritardo) mediante uncoefficiente di proporzionalita. Tale approccio conduce alla seguente equazione differenzialecon ritardo:

x′(t)x(t)

= p x(t− τ). (1.3.1)

Naturalmente la (1.3.1) puo benissimo essere modificata aggiungendo anche un termine dicontrollo logistico.

Il modello introdotto nell’esempio e senza dubbio interessante ma per certi versi limitato.Potrebbe infatti essere necessario considerare l’effetto che il numero di individui a untempo s generico (s < t) comporta sulla derivata logaritmica del numero di individui altempo t. Supponiamo che il numero di individui presenti al tempo s contribuisca al tassorelativo di crescita del numero di individui al tempo t secondo una quantita definita dauna funzione K(t, s). Appare intuitivo tener conto degli infiniti contributi che derivano daogni tempo s < t. Per fare cio si puo pensare di integrare tutti i contributi tra il tempoiniziale t0 ≤ 0 e il tempo corrente t. La relazione che ne deriva prende dunque la forma:

x′(t)x(t)

=∫ t

t0

K(t, s)x(s)ds. (1.3.2)

Tale equazione e una equazione integro-differenziale che supponiamo sia da risolversi sul-l’intervallo [0,+∞). Affinche il problema risulti ben posto e necessario assegnare un datoiniziale x(0) = x0 se t0 = 0, o una funzione iniziale x0(t) definita sull’intervallo [t0, 0] set0 < 0. Generalmente la funzione K(t, s) prende il nome di nucleo e ha il ruolo di memoriadella storia passata. Tipicamente si prediligono nuclei decrescenti, non necessariamentecontinui ma integrabili in s. Inoltre si tende a scegliere il nucleo cosicche l’integrale diventiun integrale di convoluzione ovvero a scegliere K(t, s) = K∗(t− s)5. Qualora tale modellocomprenda anche un termine di correlazione istanataneo (potenziale biologico) e un ter-mine di controllo logistico, moltiplicando ambo i membri per x(t) si ottiene la seguenteequazione:

x′(t) = k x(t)− k

εx(t)2 + x(t)

(∫ t

t0

K(t, s)g(s, x(s))ds). (1.3.3)

5Per maggiori dettagli su come l’integrale in questione possa essere visto come integrale di convoluzione,con particolare riferimento alla tipologia di nucleo trattata nel Capitolo 2 si veda il Capitolo medesimo

La (1.3.3) prende il nome di Equazione di Verhulst-Volterra e necessita di condizioni inizialiidentiche a quelle della (1.3.2).

1.3.2 Considerazioni sull’equazione di Verhulst-Volterra

L’equazione di Verhulst-Volterra e un’ equazione integro-differenziale piuttosto complicata.La casistica presenta numerosissime varianti che dipendono oltre che dai coefficienti deitermini di primo e secondo ordine anche dal numero infinito di nuclei che si possonoscegliere. La teoria presenta pochi risultati generali, ma uno particolarmente interessantee contenuto nel Teorema 1.3.1 (Teorema 4.1.2 di [BVIDE]).

Teorema 1.3.1. Siano r, R e m costanti positive. Inoltre m soddisfi la seguente:∫ t

t0

|K(t, s)|ds ≤ m 0. Allora r < x(t) < R per 0 ≤ t <∞

Dimostrazione. Supponiamo che il risultato sia falso e che t1 sia il primo valore di t per ilquale l’asserto r < x(t) < R sia violato. Per essere precisi si supponga x(t1) = R. Alloradall’equazione differenziale di ha:

x′(t1)x(t1)

≤ a− bR+R

∫ t1

t0

|K(t1, s)|ds

< a− bR+Rm < 0

il che e una contraddizione in quanto si deve avere x′(t1) ≥ 0. Ora supponiamo chex(t1) = r. Allora si avrebbe x′(t1)/r ≥ a− br −Rm > 0 che e ancora una contraddizionepoiche dovrebbe essere x′(t1) ≤ 0. E questo conclude.

Tale Teorema afferma sostanzialmente che se il nucleo appartiene a L1(t0, t) ∀t ≥ t0 e||K(t, s)||L1(t0,t) e sufficientemente piccola allora la soluzione e limitata e l’estinzione nonavviene. E importante notare come il Teorema 1.3.1 sia valido sotto ipotesi abbastanza ge-nerali sul nucleo. La trattazione esaustiva dell’equazione di Verhulst-Volterra verra svoltanel Capitolo ad essa dedicato nel caso particolare di nucleo a decadimento esponenziale.

1.4 Il sistema di Lotka-Volterra

1.4.1 Introduzione

Questo modello e il primo che viene proposto nell’ambito della descrizione di sistemi a duepopolazioni. Si suppone in particolare di considerare una interazione tra popolazioni di tipopreda-predatore. Negli anni successivi al primo conflitto mondiale, nel mare Adriatico siregistro un calo di una certa specie di pesce commestibile e un aumento del suo predatore.Si trattava di capire come mai l’assenza di pesca durante la guerra avesse favorito losviluppo del predatore e sfavorito quello della preda. La questione fu posta al matematico

Il sistema di Lotka-Volterra 7

Vito Volterra nel 1924. Egli risolse il problema proponendo il celebre modello di cui stiamotrattando. Volterra indico con x(t) il numero delle prede in funzione del tempo e con y(t) ilnumero dei predatori. Sulla base delle informazioni fornite dai biologi formulo le seguentiipotesi. Suppose il tasso relativo di crescita delle prede costante in assenza di predatorima decrescente linearmente in funzione del numero dei predatori. Viceversa assunse che ipredatori decrescessero ad un tasso costante in assenza di prede ma crescessero linearmentein funzione di x. Tali considerazioni derivano dall’idea che in assenza di predatori le predesi possano riprodurre secondo un modello di Malthus mentre che la presenza dei predatorifaccia decrescere il numero di prede mediante un fattore proporzionale al numero di incontritra prede e predatori (prodotto tra il numero di prede e quello di predatori). Al contemposi suppone che in assenza di prede, e quindi di cibo, i predatori si estinguano lentamentesempre secondo un modello Malthusiano (pero con potenziale biologico negativo) mentreche in presenza di prede si riproducano in proporzione al numero di prede mangiate (sempresupposto correlato al prodotto x y). Il modello formulato da Volterra si presentava dunquenel seguente modo: x′(t) = x(t)

(a− b y(t)

)y′(t) = y(t)

(− c+ d x(t)

) . (1.4.1)

Tale sistema a cui bisogna associare le condizioni al contorno x(0) = x0 e y(0) = y0 eun sistema di equazioni differenziali accoppiate non lineare a coefficienti costanti (tuttisupposti positivi) che prende il nome di sistema di Lotka-Volterra. Il nome di Lotka vienegeneralmente associato a quello di Volterra in quanto anche egli studio il medesimo sistemanello stesso periodo (Lotka giunse tuttavia alla formulazione in un altro ambito di ricercae del tutto indipendentemente da Volterra).

1.4.2 Analisi delle soluzioni

In questa sezione analizzeremo per quanto possibile le soluzioni del sistema (1.4.1). Ini-ziamo con l’avvisare il lettore che il sistema non ammette risoluzione in forma chiusa.Tuttavia sono possibili varie considerazioni. Per prima cosa e opportuno individuare eanalizzare gli equilibri di interesse (quelli a coordinate non negative). Il sistema presentadue possibili equilibri:

xeq1 = 0, yeq1 = 0; (1.4.2a)

xeq2 =c

d, yeq2 =

a

b. (1.4.2b)

Notiamo che entrambi gli equilibri rientrano nella zona di interesse modellistico. Lineariz-zando il sistema e valutando la matrice Jacobiana ad esso associata nei punti di equilibriosi trova facilmente che l’equilibrio (1.4.2a) e sempre instabile (autovalori opposti in segno)mentre non e possibile fare conclusioni dirette sull’equilibrio (1.4.2b) (autovalori complessie coniugati con parte reale nulla). Fortunatamente per il sistema (1.4.1) esiste un integraleprimo. Tale integrale primo e:

E(x, y) = −c ln(x) + d x− a ln(y) + b y. (1.4.3)

Si deduce dunque che gli insiemi di livello di E(x, y) sono orbite del sistema. Poiche∂2E∂x2 = c

x2 , ∂2E∂x ∂y = 0 e ∂2E

∂y2 = ay2 si ha che E e strettamente convessa in R2

+. Pertanto


Figura 1.3: Spazio delle fasi del sistema (1.4.1). Sono riportate in particolare le tre orbite corri-spondenti ai dati iniziali x0 = 2, y0 = 2 (rosso); x0 = 2.5, y0 = 2.5 (verde); x0 = 3, y0 = 3 (blu).Gli altri parametri caratterizzanti questa specifica realizzazione del sistema (1.4.1) sono: a = 1,b = 1, c = 1 e d = 1. In ordinata viene rappresentato il numero dei predatori mentre in ascissa ilnumero di prede.

l’equilibrio (1.4.2b) e un punto di minimo per E. Notando inoltre che limx→+∞E(x, y) =limy→+∞E(x, y) = limx→0+ E(x, y) = limy→0+ E(x, y) = +∞ si deduce che gli insiemi dilivello di E sono curve chiuse regolari che contengono il punto di equilibrio al loro interno.Cio comporta che:

Microteorema 1.4.1. Qualora i dati iniziali siano strettamente positivi la soluzione delsistema di Lotka-Volterra e periodica.

Come gia preannunciato per tali soluzioni non esiste una espressione esplicita. Il lettore etuttavia invitato a vedere gli andamenti di alcune orbite nel piano delle fasi (Figura 1.3)e di una soluzione in funzione del tempo (Figura 1.4) nel caso in cui i parametri sianoa = b = c = d = 1. Per quanto possa sembrare strano, visto che il sistema e e continuaad essere oggetto di studio da parte di grandi esponenti della comunita matematica, none tuttora nota alcuna formula esplicita che esprima il periodo del sistema in funzionedei parametri e del dato iniziali. Vi sono formule integrali che coinvolgono la soluzionenumerica ma cio esula dagli scopi del nostro lavoro. Vorremmo tuttavia mettere in risaltocome sia possibile determinare il valor medio delle componenti della soluzioni nell’arco di


Figura 1.4: Soluzione numerica del sistema (1.4.1). In tale grafico vengono illustrati gli andamentidel numero di prede (blu) e di predatori (rosso) in funzione del tempo. I parametri caratterizzantiil modello in questo caso specifico sono a = 1, b = 1, c = 1 e d = 1 mentre la coppia di dati inizialie x0 = y0 = 2.5. In ordinata si trovano i valori numerici delle due variabili mentre in ascissa iltempo.

un periodo. In particolare, detto T il periodo, vale il seguente Teorema (Cfr. [PS2, pagg.291-293]):

Teorema 1.4.2.1T

∫ T

0y(t)dt =

a

b= yeq2 ; (1.4.4a)

1T

∫ T

0x(t)dt =

c

d= xeq2 . (1.4.4b)

Durante il nostro lavoro siamo riusciti ad ottenere alcuni risultati relativi al sistema diLotka-Volterra. Tali risultati sono raccolti nella sezione 1.4.3.

1.4.3 Teoremi sul sistema di Lotka-Volterra

Si consideri il sistema di Lotka-Volterra:

Σ :

x′ = a x− b x y

y′ = −c y + d x y(1.4.5)

Supporremo senza perdita di generalita all’interno del quadro modellistico che i parametria, b, c, d e i dati iniziali x0, y0 siano tutti strettamente positivi. Il primo Teorema cheproponiamo riguarda il legame delle componenti x e y delle soluzioni al variare di alcuniparametri e sotto opportune condizioni:

Teorema 1.4.3. Si consideri il sistema (1.4.5) dipendente dai parametri a, b, c, d. Allo-ra, fissati i parametri a, b e c, il dato iniziale y0 e una costante arbitraria k positiva, lacomponente y(t) della soluzione e la medesima per tutti i sistemi tali che valga d x0 = k.Inoltre la componente x(t) della soluzione dipende dalla particolare scelta di d e di x0 soloin termini di un coefficiente moltiplicativo.


Dimostrazione. Si consideri il sistema (1.4.5), si ricavi x(t) dalla seconda equazione e sisostituisca nella prima dopo averne calcolato la derivata.

x =y′ + c y

dy

x′ =d y y′′ − d(y′)2

(d y)2

d y y′′ − d(y′)2

(d y)2= a

y′ + c y

d y− b

y′ + c y

d yy

Ricordandosi che si e imposto d > 0, e x0, y0 > 0 e considerato che nel sistema di Lotka-Volterra le soluzioni sono periodiche, e se partono nel primo quadrante non si annullanomai, si puo semplificare per d y e svolgendo i conti si perviene all’equazione differenziale:

y′′ =(y′)2

y+ (a− b y)(c y + y′)

a cui vanno associate le condizioni iniziali y(0) = y0 e y′(0) = −c y0 + d x0 y0. Taleequazione si puo riscrivere come un sistema di equazioni differenziali ordinarie non linearedel primo ordine del tipo w′ = f(w) con f localmente lipshitziana in w uniformemente int. Pertanto per il Teorema di Cauchy-Lipshitz6 esiste una e una sola soluzione y(t). Poichel’equazione non dipende dal parametro d se non nell’assegnamento di y′(0) risulta evidenteche scegliendo i sistemi per cui d x0 e costante i dati iniziali non variano e di conseguenzala componente x(t) del sistema e sempre la stessa. Per la seconda parte si ricavi y(t) nellaprima equazione e la si sostituisca nella seconda dopo averne calcolato la derivata.

y =a x− x′

b x

y′ =b(x′)2 − b x x′′

(b x)2

b(x′)2 − b x x′′

(b x)2= −ca x− x′

b x+ d x

ax− x′

b x

e proseguendo come in precedenza si perviene all’equazione differenziale:

x′′ =(x′)2

x+ (c− d x)(a x− x′)

Si supponga ora di considerare due sistemi Σ1 e Σ2:

Σ1 :

x′1 = a x1 − b x1 y1

y′1 = −c y1 + d1 x1 y1

Σ2 :

x′2 = a x2 − b x2 y2

y′2 = −c y2 + d2 x2 y2

6vedasi appendice.


Si scelga y1(0) = y2(0) = y0, e x1(0) d1 = x2(0) d2. Per quanto detto in precendeza si hache y1(t) ≡ y2(t) ∀t > 0. Mostriamo ora che x1(t) = d2

d1x2(t). L’equazione differenziale e

le condizioni iniziali per x2(t) sono:

x′′2 −(x′2)

2

x2+ (d2x2 − c)(a x2 − x′2) = 0

x2(0) = x2

x′2(0) = a x2 − b x2y0

Supponiamo che x2(t) = d1d2x1(t). Sostituendo nell’equazione differenziale e nelle condi-

zioni al contorno si ha:

x′′1 −(x′1)

2

x1+ (d1x1 − c)(a x1 − x′1) = 0

x1(0) = x1

x′1(0) = a x1 − b x1y0

Che e proprio l’equazione differenziale per x1(t).

Risultati numerici utili a comprendere il significato del Teorema 1.4.3 sono riportati inFigura 1.5 e 1.6. In particolare vengono proposte le orbite di due sistemi e le corrispondentisoluzioni in funzione del tempo. I parametri mantenuti fissi sono: a = b = c = 1 e il datoiniziale y0 = 2. Si e scelto come costante arbitraria k > 0 il valore 2. Il sistema 1 presentacome altri parametri x01 = d1 =

√2 mentre il sistema 2 ha x02 = 2, d2 = 1. In base a

questi valori e a quanto affermato nel Teorema il lettore puo trovare riscontro visivo delfatto che y1(t) = y2(t) ∀t e che x2(t) = d1

d2x1(t) =

√2x1(t). Il Teorema che segue puo

considerarsi il duale del precedente.

Teorema 1.4.4. Si consideri il sistema (1.4.5) dipendente dai parametri a, b, c, d. Allora,fissati i parametri a, c e d, il dato iniziale x0 e una costante arbitraria k positiva, lacomponente x(t) della soluzione e la medesima per tutti i sistemi tali che valga b y0 = k.Inoltre la componente y(t) della soluzione dipende dalla particolare scelta di b e di y0 soloin termini di un coefficiente moltiplicativo.

Dimostrazione. La dimostrazione segue lo stesso schema di quella del Teorema (1.4.3) eviene pertanto omessa.

Corollario 1.4.5. Si consideri il sistema (1.4.5) dipendente dai parametri a, b, c, d. Al-lora fissati i parametri a e c e due costanti arbitrarie k1 e k2 positive, tutte le soluzionicorrispondenti a sistemi i cui dati iniziali soddisfino b y0 = k1 e d x0 = k2 hanno lo stessoperiodo T . Inoltre le componenti x(t) e y(t) della soluzione variano al variare delle sceltedi b, y0, d, x0 solo mediante due costanti moltiplicative.

Dimostrazione. La dimostrazione e banale conseguenza dei Teoremi (1.4.3) e (1.4.4).

Il prossimo Teorema riguarda il particolarissimo caso in cui i parametri del sistema e i datiiniziali siano uguali a coppie.


Figura 1.5: Simulazione numerica di due sistemi atta a comprendere meglio l’enunciato del Teo-rema 1.4.3. In questa figura vengono proposte le orbite di due sistemi con identici parametri a = 1,b = 1, c = 1 e con identico dato iniziale y0 = 2. La costante arbitraria scelta e k = 2. In rossole orbite del sistema con x01 = d1 =

√2 mentre in blu quelle del sistema con x02 = 2, d2 = 1. In

ordinata viene rappresentato il numero di predatori mentre in ascissa il numero di prede.

Figura 1.6: Simulazione numerica di due sistemi atta a comprendere meglio l’enunciato del Teo-rema 1.4.3. In questa figura vengono proposte le soluzioni in funzione del tempo di due sistemi conidentici parametri a = 1, b = 1, c = 1 e con identico dato iniziale y0 = 2. La costante arbitrariascelta e k = 2. A sinistra le orbite del sistema con x01 = d1 =

√2 mentre a destra quelle del sistema

con x02 = 2, d2 = 1. Nelle figure sono tracciati i grafici in funzione del tempo delle prede (blu)e dei predatori (rosso). Si noti in blu tratteggiato una opportuna dilatazione della variabile x(t)relativa al sistema con x01 = d1 =

√2 che coincide esattamente con la variabile x(t) del sistema con

x02 = 2, d2 = 1. Si noti inoltre che come enunciato dal teorema 1.4.3 la soluzione per la variabiley(t) e la medesima per entrambi i sistemi. In ordinata viene rappresentato il numero di prede epredatori mentre in ascissa il tempo.

Teorema 1.4.6. Si consideri il sistema (1.4.5) dipendente dai parametri a, b, c, d. Allorase a = c , b = d e x0 = y0 la soluzione ha la forma:(

x(t)y(t)

)=(

x(t)x(−t)

)=(y(−t)y(t)

)Dimostrazione. Si consideri l’equazione differenziale associata alla y ricavata nella dimo-strazione del Teorema (1.4.3). Si supponga y(t) = x(−t). L’equazione dopo aver sostituitola x in luogo della y diventa:

x′′ =(x′)2

x+ (a− b x)(c x− x′) =

(x′)2

x+ (c− d x)(a x− x′)

a cui vanno associate le condizioni al contorno x0 = x(0) = y(0) = y0 e x′(0) = −y′(0) =c y0 − d x0 y0 = −a y0 + b x0 y0 ovvero l’equazione differenziale generica per la x e le con-dizioni al contorno ad essa associate ricavate nella dimostrazione del Teorema (1.4.3).

Il Teorema che segue puo a tutti gli effetti considerarsi una estensione del Teorema 1.4.2.

Teorema 1.4.7. Si consideri il sistema (1.4.5) dipendente dai parametri a, b, c, d. Allorase T e il periodo e k ∈ N\0 allora vale:

limt→∞

1t

∫ t

0y(s) ds =

1kT

∫ kT

0y(s)ds =

a

b

Dimostrazione. Dividendo la prima equazione del sistema per x e integrando tra 0 e t siottiene:

ln(x(t)x0

)= a t− b

∫ t

0y(s)ds

Allora dividendo per t e riarrangiando i termini si ha:

1t

∫ t

0y(s) ds =

a

b+ ln

((x(t)x0

)− 1b t

)

da cui segue immediatamente il risultato per t = kT . Tenendo inoltre conto che le orbitedel sistema di Lotka-Volterra (1.4.5) sono periodiche e sempre nel primo quadrante perogni coppia di dati iniziali strettamente positivi esistono certamente due valori xmin exmax tali da soddisfare 0 < xmin ≤ x(t) ≤ xmax <∞ ∀t ≥ 0. Dunque considerato che illogaritmo e le potenze sono funzioni crescenti e lecito scrivere:

a

b+ ln

((xmin

x0

)− 1b t

)≤ 1t

∫ t

0y(s) ds ≤ a

b+ ln

((xmax

x0

)− 1b t

)

e passando al limite si conclude.

Il Teorema che segue e duale rispetto al Teorema 1.4.7



limt→∞

1t

∫ t

0x(s) ds =

1kT

∫ kT

0x(s) ds =

c

d

Dimostrazione. La dimostrazione e sostanzialmente identica a quella del Teorema (1.4.7)e pertanto viene lasciata al lettore.

L’ultimo Teorema che proponiamo lega il prodotto delle medie delle variabili alla mediadel prodotto delle variabili stesse.


limt→∞

1t

∫ t

0x(s)y(s) ds =

1kT

∫ kT

0x(s)y(s) ds =

ac

bd

Dimostrazione. Siano x e y le medie integrali delle variabili su un periodo T . Allora sipuo scrivere:

1t

∫ t0 x(s)y(s) ds =

= 1t

∫ t0 (x(s)− x)(y(s)− y) ds+ x

t

∫ t0 y(s) ds+ y

t

∫ t0 x(s) ds− xy =

= 1bdt

∫ t0

y′(s)y(s) (by(s)− a) ds+ x

t

∫ t0 y(s) ds+ y

t

∫ t0 x(s) ds− xy =

= 1bdt

(by(t)− by0 − a ln

(y(t)y0

))+ x

(ab + ln

((x(t)x0

)− 1bt

))+

+y(

cd + ln

((y(t)y0

) 1d t

))− xy =

= acbd + y(t)−y0

dt + ln((

x(t)x0

)− cdbt

)

Da cui passando al limite o considerando che x(kT ) = x0 e y(kT ) = y0 si conclude.

In Figura 1.7 viene mostrato l’andamento del prodotto delle medie, della media del pro-dotto e del prodotto delle medie sul periodo (costante) in funzione del tempo nel caso incui a = 1, b = 1.5, c = 4, d = 1.5, x0 = 1, y0 = 3. Tale Figura risulta di fondamentaleimportanza per comprendere al meglio il Teorema 1.4.9. Inoltre, ripetendo le prove nu-meriche al variare dei parametri e confrontandole con i periodo ottenuti sempre per vianumerica, sembra ci sia sufficiente evidenza per poter formulare la seguente congettura.

Congettura 1.4.1. Il periodo T e il piu piccolo t > 0 per il quale si verifica:

1t

∫ t

0x(s)y(s) ds =

(1t

∫ t

0x(s)ds

)(1t

∫ t

0y(s)ds

)=ac

bd

Con riferimento alla Figura 1.7 si nota come il primo punto in cui tutti e tre i graficicoincidono sia circa 4, 3. Tale valore coincide numericamente con il valore del periodo. Sinoti inoltre che, anche se sembra che vi sia un punto antecedente in cui i grafici analizzati siintersechino (circa 0, 2) cio viene smentito da un ingrandimento della Figura. Per mostrarecon maggior forza tale concetto si vedano Figura 1.8.

Il sistema di Lotka-Volterra con controllo logistico 15

Figura 1.7: Simulazione numerica atta a supportare la congettura 1.4.1. Sono riportati gli anda-menti in funzione del tempo del prodotto delle medie integrali (blu) delle variabili x(t) e y(t), dellamedia integrale del prodotto delle variabili (rosso) e di a c

b d (verde). I valori numerici dei parametrie dei dati iniziali relativi a questa simulazione sono: a = 1, b = 1.5, c = 4, d = 1.5, x0 = 1 e y0 = 3.In ordinata viene riportato il valore numerico delle tre funzioni citate mentre in ascissa il tempo.

Figura 1.8: Ingrandimenti relativi ai grafici di Figura 1.7 atti a chiarire meglio come il primovalore in cui dovrebbe sussistere la relazione della congettura 1.4.1 sia esattamente il periodo. Ildubbio che potesse valere la relazione in corrispondenza di un valore temporale minore del periodo,che la Figura 1.7 poteva suscitare, viene chiaramente fugato.


1.5 Il sistema di Lotka-Volterra con controllo logistico

1.5.1 Introduzione

Nella sezione 1.4 abbiamo introdotto il sistema di Lotka-Volterra come modello per descri-vere un’interazione di tipo preda-predatore tra due popolazioni. Naturalmente, come giaosservato nel caso di singola popolazione tale modello risulta suscettibile di raffinamenti.Il raffinamento proposto in questa parte prevede di tener conto della possibile limitatez-za di risorse dell’habitat aggiungendo dei termini di controllo logistico. Se il lettore nonavesse letto la sezione 1.2 e invitato a farlo per comprendere cosa si intenda per controllologistico e come tale controllo venga modellizzato. Modificando il sistema (1.4.1) mediantei parametri di controllo logistico sopra citati si perviene al sistema (1.5.1) seguente: x′(t) = x(t)

(a− b y(t)− g x(t)

)y′(t) = y(t)

(− c+ d x(t)− γ y(t)

) (1.5.1)

Tale sistema e un sistema di equazioni differenziali del primo ordine non lineari a coefficienticostanti. Al fine di garantire il significato modellistico si suppone che tutti i parametri sianostrettamente positivi o al limite che uno dei due controlli logistici sia nullo. Naturalmenteaffinche il problema risulti ben posto e necessario assegnare una coppia di dati iniziali x(0)e y(0), anch’essi supposti non negativi ai fini modellistici. Infine notiamo che anche inquesto caso soluzioni che, come richiesto, abbiano dati iniziali non negativi si mantengononon negative.

1.5.2 Analisi delle soluzioni

Come gia per il sistema di Lotka-Volterra senza controllo logistico, anche in questo casonon e possibile ottenere una soluzione in forma chiusa. Il sistema in questione presenta4 punti di equilibrio di cui due con entrambe le coordinate sempre non negative, unocon entrambe le coordinate non negative in dipendenza dai parametri e l’ultimo con unacoordinata sempre negativa. In questa breve trattazione ci focalizzeremo come piu voltesottolineato solo sugli equilibri di interesse modellistico ovvero quelli a coordinate nonnegative. I tre equilibri in questione sono:

xeq1 = 0, yeq1 = 0 (1.5.2a)

xeq2 =a

g, yeq2 = 0 (1.5.2b)

xeq3 =b c+ a γ

b d+ γ g, yeq3 =

a d− c g

b d+ γ g(1.5.2c)

In particolare gli equilibri (1.5.2a) (estinzione di entrambe le specie) e (1.5.2b) (estinzionedei soli predatori) risultano sempre interessanti dal punto di vista modellistico mentrel’equilibrio (1.5.2c) (coesistenza di entrambe le specie) assume interesse solo se a d ≥ g c.Si noti che nel caso in cui valga esattamente a d = g c gli equilibri (1.5.2b) e (1.5.2c)coincidono. Se infatti si chiama xeq2 = ϑ si puo scrivere a = ϑ g. L’uguaglianza a d = g cimplica che yeq3 = 0 = yeq2 . Inoltre dalle due uguaglianze imposte si ricava agevolmente che

Il sistema di Lotka-Volterra con memoria 17

c = d ϑ. Infine sfruttando le relazioni introdotte si ottiene facilmente che anche xeq3 = ϑe pertanto i due equilibri coincidono.La linearizzazione del sistema porta alla seguente matrice Jacobiana:

J =[a− b y − 2g x −b x

d y −c+ d x− 2γ y

](1.5.3)

Valutando la matrice (1.5.3) nelle coordinate relative all’equilibrio (1.5.2a) si nota subitocome tale matrice abbia due autovalori di segno opposto e come pertanto l’equilibrioconsiderato sia sempre instabile. Si puo dunque enunciare il seguente:

Microteorema 1.5.1. Se ambedue le specie hanno un numero iniziale di individui stret-tamente positivo non e possibile che si estinguano entrambe.

Valutando invece la matrice (1.5.3) nelle coordinate dell’equilibrio (1.5.2b) si scopre che ta-le equilibrio e localmente asintoticamente stabile se d a < g c mentre instabile se d a > g c.Per l’uguaglianza e necessaria una analisi piu approfondita. Notiamo tuttavia l’equilibrio(1.5.2b) risulta localmente asintoticamente stabile finche l’equilibrio (1.5.2c) presenta unacoordinata negativa, ovvero non ha significato dal punto di vista modellistico. Come delresto ci si potrebbe aspettare, valutando la matrice (1.5.3) nelle coordinate dell’equilibrio(1.5.2c) si vede che se tale equilibrio e a coordinate non negative risulta localmente asinto-ticamente stabile7. Invitiamo ancora il lettore a prestare attenzione al fatto che l’equilibrio(1.5.2c) coincide con l’equilibrio (1.5.2b) qualora abbia la coordinata yeq3 = 0 ovvero per ilprimo valore di yeq3 che lo rende significativo dal punto di vista modellistico. Sintetizziamoquanto detto nei seguenti riultati:

Microteorema 1.5.2. Se l’equilibrio (1.5.2c) (coesistenza di entrambe le specie) ha si-gnificato dal punto di vista modellistico allora e localmente asintoticamente stabile mentrel’equilibrio (1.5.2b) (estinzione dei soli predatori) e instabile. Viceversa se l’equilibrio(1.5.2c) non ha rilevanza modellistica allora l’equilibrio (1.5.2b) e localmente asintotica-mente stabile.

Come gia fatto nella sezione 1.4 per il sistema di Lotka-Volterra semplice vengono propostialcuni grafici che illustrano la soluzione del sistema di Lotka-Volterra con controllo logisticoin corrispondenza a scelte particolari dei parametri e dei dati. Nelle Figure 1.9 e 1.10 sonoriportati rispettivamente la traiettoria nello spazio delle fasi e la soluzione in funzione deltempo qualora si tenda asintoticamente verso l’equilibrio (1.5.2b). Nelle Figure 1.11 e 1.12sono messi in evidenza gli stessi aspetti ma relativi al caso in cui l’equilibrio asintoticamentestabile sia quello della (1.5.2c).

1.6 Il sistema di Lotka-Volterra con memoria

1.6.1 Introduzione

Per introdurre il sistema di Lotka-Volterra con memoria si potrebbero fare delle conside-razioni analoghe a quelle fatte nella sezione 1.3 in merito all’introduzione dell’equazione

7Per la dimostrazione si e sfruttato il Teorema A.1.1


Figura 1.9: Spazio delle fasi del sistema (1.5.1) qualora l’equilibrio (1.5.2b) risulti attrattivo.Viene riportata l’orbita corrispondente ai dati iniziali x0 = 2.5, y0 = 2.5 (rosso). Gli altri parametricaratterizzanti questa specifica realizzazione del sistema (1.5.1) sono: a = 1, b = 1, c = 1, d = 1,g = 2 e γ = 1. Si noti come i predatori tendano ad estinguersi mentre il numero di prede tendaa stabilizzarsi su una quantita costante. In ordinata viene rappresentato il numero dei predatorimentre in ascissa il numero di prede.

Figura 1.10: Soluzioni in funzione del tempo del sistema (1.5.1) qualora l’equilibrio (1.5.2b) risultiattrattivo. Vengono riportati i grafici in funzione del tempo del numero di prede (blu) e predatori(rosso) corrispondenti ai dati iniziali x0 = 2.5, y0 = 2.5. Gli altri parametri caratterizzanti questaspecifica realizzazione del sistema (1.5.1) sono: a = 1, b = 1, c = 1, d = 1, g = 2 e γ = 1. Sinoti come i predatori tendano ad estinguersi mentre il numero di prede tenda a stabilizzarsi suuna quantita costante. In ordinata viene rappresentato il numero dei predatori mentre in ascissail numero di prede.


Figura 1.11: Spazio delle fasi del sistema (1.5.1) qualora l’equilibrio (1.5.2c) risulti attrattivo.Viene riportata l’orbita corrispondente ai dati iniziali x0 = 2.5, y0 = 2.5 (rosso). Gli altri parametricaratterizzanti questa specifica realizzazione del sistema (1.5.1) sono: a = 1, b = 1, c = 1, d = 5,g = 1 e γ = 1. Si noti come sia la quantita di predatori che di prede tendano a stabilizzarsi suun valore costante. In ordinata viene rappresentato il numero dei predatori mentre in ascissa ilnumero di prede.

Figura 1.12: Soluzioni in funzione del tempo del sistema (1.5.1) qualora l’equilibrio (1.5.2c) risultiattrattivo. Vengono riportati i grafici in funzione del tempo del numero di prede (blu) e predatori(rosso) corrispondenti ai dati iniziali x0 = 2.5, y0 = 2.5 (rosso). Gli altri parametri caratterizzantiquesta specifica realizzazione del sistema (1.5.1) sono: a = 1, b = 1, c = 1, d = 5, g = 1 e γ = 1.Si noti come sia la quantita di predatori che di prede tendano a stabilizzarsi su un valore costante.In ordinata viene rappresentato il numero dei predatori mentre in ascissa il numero di prede.


integro-differenziale scalare di Verhulst-Volterra. Tuttavia sembra piu istruttivo seguireil ragionamento proposto dallo stesso Volterra nella formulazione originale del sistema.Si supponga di dover modellizzare un’ interazione tra specie di tipo preda-predatore. Siassuma inoltre che, perlomeno per quanto riguarda i predatori, la distribuzione per eta delnumero di individui rimanga costante nel tempo e che φ(ξ)dξ sia il rapporto tra il numerodi individui di eta compresa tra ξ e ξ + dξ e il numero totale di predatori N2(t). Allora ilnumero di predatori che al tempo t ha almeno ϑ anni si puo ottenere tramite la seguenteespressione:

N2(t)∫ ∞

ϑφ(ξ)dξ = N2(t) f(ϑ), f ≥ 0 , f 6≡ 0

Tra gli N2(t) individui esistenti al tempo t, ce ne saranno N2 f che esistevano gia unprecedente tempo τ . Chiamato N1(t) il numero di prede in funzione del tempo, dalla leggedi azione di massa si puo vedere che la quantita di nutrimento assorbito nell’intervallo(τ, τ + dτ) dai predatori che esistevano ad entrambi i tempi τ e t e:

γ f(t− τ)N1(τ)N2(t)dτ

ove γ e una opportuna costante di proporzionalita. Si puo inoltre pensare che talenutrimento produca un incremento di predatori nell’intervallo (t, t+ dt) pari a:

ψ(t− τ)dtγ f(t− τ)N1(τ)N2(t)dτ, ψ ≥ 0, ψ 6≡ 0

Allora, sommando i contributi infinitesimi di tali incrementi tra −∞ e t la secondaequazione del sistema di Lotka-Volterra (1.4.1) si trasforma nella seguente:

dN2 = −ε2N2dt+N2(t)dt∫ t

−∞ψ(t− τ)γ f(t− τ)N1(τ)dτ

ove ε2 svolge il ruolo del parametro c del sistema (1.4.1). Volterra scrisse dunque il sistema:N ′

1(t) = N1(t)(ε1 − γ1N2(t)

)N ′

2(t) = −ε2N2 +N2(t)∫ t−∞ ψ(t− τ)γ f(t− τ)N1(τ)dτ

ove ε1 e γ1 hanno rispettivamente il ruolo dei parametri a e b contenuti nel sistema (1.4.1).La formulazione originale di Volterra puo essere ulteriormente raffinata aggiungendo itermini di controllo logistico e un termine di interazione con memoria a primo membro.In particolare siamo interessati alla seguente formulazione che prende il nome di sistemaintegro-differenziale di Lotka-Volterra: x′(t) = x(t)

(a− bx(t)− gy(t)−

∫ tt0K1(t, s)y(s)ds

)y′(t) = y(t)

(− α− βy(t) + γx(t) +

∫ tt0K2(t, s)x(s)ds

) (1.6.1)

Nella (1.6.1) i nomi delle funzioni descriventi il numero di prede e predatori in funzionedel tempo sono diversi da quelli utilizzati da Volterra nella formulazione originale e ripro-posti in precedenza. Tale scelta e stata fatta per seguire una nomenclatura piu consona inrelazione ai modelli 1.4 e 1.5 nonche a [BVIDE]. Inoltre sono stati cambiati anche i nomi

dei parametri del sistema. In tale scelta si e seguita una discontinuita rispetto ai modelliprecedenti e a [BVIDE] per evidenziare il ruolo duale dei parametri nelle due equazioni. Inparticolare nella prima equazione, quella che descrive l’evoluzione della variabile x(t) ovve-ro del numero di prede, vengono utilizzate lettere romane mentre nella seconda equazione,relativa all’evoluzione del numero dei predatori, lettere greche. I parametri a e α hannoora il ruolo di potenziale biologico, i parametri b e β di controllo logistico, i parametri g e γdi interazione istantanea preda-predatore e infine al nucleo proposto da Volterra vengonosostituiti due nuclei generici non negativi (uno per ogni equazione) indicati con K1(t, s) econ K2(t, s).

1.6.2 Considerazioni sul sistema di Lotka-Volterra con memoria

Lo studio del sistema (1.6.1) sara l’oggetto principale del Capitolo 4. Tale studio verrasvolto nel caso particolare in cui i nuclei K1 e K2 siano a decadimento esponenziale comequelli presentati nel Capitolo 2 e gli integrali siano di convoluzione. La teoria presentapochissimi risultati in merito e sarebbe interessante trarre degli andamenti anche soloqualitativi circa la possibile esistenza di cicli limite, orbite periodiche o punti di equilibrio.Un risultato sul sistema (1.6.1) proposto da Burton in [BVIDE] come Teorema 4.1.3 vienequi riproposto nel Teorema seguente.

Teorema 1.6.1. Siano Ki(t, s) ≥ 0 e continui per t0 ≤ s ≤ t < ∞. Sia β = 0 e tutti glialtri parametri costanti positive. Se esiste ε > 0 che soddisfa:

α >a

b

∫ t

t0

K2(t, s)ds+ ε

per ogni t ≥ t0. Allora tutte le soluzioni contenute nel primo quadrante e che soddisfanox(t) ≤ a

b sull’intervallo iniziale [t0, 0] sono limitate.

Dimostrazione. Per prima cosa si noti che x′ < 0 se x > ab , quindi, x(t) ≤ a

b se t > t0. Poisi consideri la seguente catena di disuguaglianze:

γ x′ + g y′ ≤ γ x (a− b x− g y) + g y (−α+ γ x+∫ t

t0

K2(t, s)x(s)ds)

≤ γ x (a− b x) + g y (−α+a

b

∫ t

t0

K2(t, s)ds)

≤ γ x(a− b x)− g ε y

L’ultimo termine e di conseguenza anche i precedenti sono negativi se:

g ε y > a γ x− b γ x2

Quindi se la curva γ x + g y = costante si trova sopra la parabola g ε y = a γ x − b γ x2

allora (γ x+ g y)′ ≤ 0 da cui segue che anche la y e limitata e questo conclude.

Capitolo 2

Nuclei a decadimento esponenziale

Oggetto principale di questo Capitolo e l’introduzione della famiglia di nuclei che verra uti-lizzata nei Capitoli 3 e 4 in merito all’ equazione integro-differenziale di Verhulst-Volterra eal sistema integro-differenziale di Lotka-Volterra con memoria. Un’ apposita sezione vienededicata a come questa tipologia di nucleo consenta di passare da un sistema integro-differenziale di n equazioni in n funzioni incognite ad un sistema differenziale autonomodi 2n equazioni in 2n incognite. Viene inoltre fatto un breve cenno alla convergenza nelsenso delle distribuzioni di tale nucleo alla delta di Dirac, sottolineandone il significatomodellistico. Infine vengono proposti alcuni esempi di convoluzione tra il nucleo e funzionielementari.

2.1 Introduzione

Nell’equazione integro-differenziale di Verhulst-Volterra e nel sistema di Lotka-Volterracon memoria, i termini integrali hanno la seguente forma:∫ t

t0

z(t, s)w(s)ds (2.1.1)

Mostriamo come tali integrali, sotto opportune ipotesi, possano essere visti come inte-grali di convoluzione. Tale interpretazione e importante in quanto consente di sfruttare,nel seguito della trattazione, i risultati gia noti su tale operazione. Ricordiamo che perconvoluzione in R tra due funzioni f e g si intende la seguente:

(f ∗ g)(t) =∫

Rf(t− s)g(s)ds (2.1.2)

Gli integrali della forma (2.1.1), in genere non possono essere scritti nella forma (2.1.2).Tuttavia, limitando l’attenzione ai soli nuclei z(t, s) della forma z(t− s) e supponendo chetali nuclei siano non nulli sul solo semiasse non negativo, si riesce a riportare la (2.1.1) allaforma (2.1.2). Infatti, introducendo delle opportune funzioni di Heaviside si puo scrivere:∫ t

t0

z(t− s)w(s)ds =∫ +∞

−∞z(t− s)H(t− s)w(s)H(s− t0)ds = ((z H) ∗ (wHt0)) (2.1.3)

23

24 Nuclei a decadimento esponenziale

Ove Ht0 indica la traslata della funzione di Heaviside che abbia il punto di discontinuitain t0. In particolare, dato che negli integrali delle equazioni o dei sistemi di interesseuna delle due funzioni e quella incognita, il cui dominio di definizione e inferiormentelimitato da t0, si puo pensare di identificare il nucleo con una funzione della forma z Hmentre la funzione incognita x con xHt0

1. Poiche la convoluzione e ben definita quandouna funzione “bilancia” l’altra, e ragionevole supporre che tale operazione sia ben definitaogniqualvolta la funzione incognita non cresca troppo rapidamente2. Una dimostrazionerigorosa della buona posizione del problema di convoluzione nel caso trattato e immediataqualora t0 = 0 e x ∈ L1

loc(R) grazie al Teorema A.2.3. Sempre grazie a tale Teorema sideduce la buona posizione del problema qualora ci si limiti a tempi t ≤ T < +∞. Infattibasta considerare la limitatezza del supporto nella definizione della funzione incognita (sipassa da x(t) a x(t)H(t − t0)H(T − t)) e sfruttare il fatto che a quel punto la funzioneincognita ha supporto compatto. Sempre mediante considerazioni analoghe, ma grazie alTeorema A.2.4 si deduce che il supporto della convoluzione e nel primo caso [0,+∞) e nelsecondo [t0, T ].

2.2 Il nucleo

Il nucleo ha il ruolo di determinare la storia della variabile pesando in modo diverso, tra-mite la convoluzione, i valori della stessa. In genere si e interessati a nuclei non negativi edecrescenti sul semiasse destro. Una scelta di questo tipo ha il significato di tenere maggior-mente in considerazione valori della variabile recenti nel tempo mentre di “dimenticare”valori della variabile lontani nel tempo. Naturalmente tanto piu e rapido il decadimentodel nucleo, tanto minore e la memoria dei valori passati della variabile. Il caso in cui lamassa del nucleo rimanga costante e si concentri vicino all’origine corrisponde, a livello dimodellizzazione, al tenere in considerazione solo i valori piu recenti della variabile stessa.Tale procedimento, portato al limite, comporta che la storia della variabile si riduca alsolo valore corrente ovvero che vi sia assenza di memoria. Dal punto di vista matematico,tale proprieta, si formalizza con la verifica che il nucleo converga, al variare di un qualcheparametro, alla delta di Dirac nel senso delle distribuzioni.Il particolare tipo di nucleo che si e scelto di utilizzare nel seguito della trattazione e che vie-ne dettagliatamente presentato in questo Capitolo, in accordo con quanto precedentementescritto, ha la seguente forma:

K(ξ) =k

rexp

− ξ

r

H(ξ) (2.2.1)

ove H(ξ) e la cosiddetta funzione di Heaviside. Il nucleo prototipale di tale famiglia(k = r = 1) viene mostrato in Figura 2.1 Come si puo vedere tale tipologia di nucleodipende da due parametri: k e r. Il ruolo di questi parametri e il loro nome viene descrittoqui di seguito:

• r: parametro di dispersione del nucleo (PDN). E il parametro che determina lavelocita con cui il nucleo (decrescente sul semiasse positivo delle ascisse) tende a

1Naturalmente se la funzione incognita e y si identifichera con y Ht02Tale supposizione e peraltro ben supportata dal fatto che i casi di interesse modellistico sono quelli per

i quali la popolazione non diverge

Da sistemi integro-differenziali a sistemi autonomi 25

Figura 2.1: Grafico del nucleo (2.2.1) nel caso base in cui k = r = 1. In ordinata il valore numericodi tale nucleo mentre in ascissa la variabile di dipendenza.

zero. Nella trattazione viene chiamato di dispersione perche modifica il modo in cuila massa (area sottesa alla curva) si distribuisce lungo l’asse medesimo. Alcuni esempinel caso in cui k = 1 vengono riportati in Figura 2.3. Si noti che tale parametro noninfluenza la norma L1(R) del nucleo.

• k: parametro di massa del nucleo (PMN). Questo parametro e quello che determinala norma integrale L1(R) del nucleo. Infatti coincide esattamente con l’integrale sututto R del nucleo medesimo. Geometricamente risulta un parametro di dilatazionenella direzione delle ordinate. Alcuni semplici esempi nel caso in cui r = 1 sonoriportati in Figura 2.2.

2.3 Da sistemi integro-differenziali a sistemi autonomi

In questa sezione viene spiegato il motivo principale per il quale si e scelta, all’internodella trattazione, la particolare famiglia di nuclei (2.2.1). Infatti per tale famiglia dinuclei (2.2.1) e possibile passare da una equazione integro-differenziale o da un sistemaintegro-differenziale, in cui l’integrale sia una convoluzione3 tra una funzione dell’incognitae il nucleo stesso, ad un sistema differenziale autonomo4. Per prima cosa si consideri

3A rigore, come gia sottolineato in 2.5 gli integrali a cui siamo interessati non sono convoluzioni veree proprie ma si possono ricondurre a convoluzioni introducendo delle funzioni di Heaviside. Per maggioridettagli si veda 2.5.

4In realta non e sempre detto che tale sistema sia autonomo. E sempre possibile passare ad un sistemapiu semplice che sia privo della convoluzione. L’ipotesi che garantisce il fatto che il sistema sia ancheautonomo viene formulata nel seguito di questa sezione. Nei casi trattati in questo elaborato i sistemisaranno comunque sempre autonomi.


Figura 2.2: Grafico atto a chiarire meglio come cambi il nucleo (2.2.1) al variare del parametrok chiamato anche PMN. I grafici sono relativi a diversi valori di k e ugual PDN. In particolare siha: r = 1 e k = 1 (blu), k = 2 (verde), k = 0.5 (rosso). In ordinata si trovano i valori numerici delnucleo mentre in ascissa la variabile di dipendenza.

Figura 2.3: Grafico atto a chiarire meglio come cambi il nucleo (2.2.1) al variare del parametror chiamato anche PDN. I grafici sono relativi a diversi valori di r e ugual PMN. In particolare siha: k = 1 e r = 1 (blu), r = 2 (verde), r = 0.5 (rosso). In ordinata si trovano i valori numerici delnucleo mentre in ascissa la variabile di dipendenza.

Da sistemi integro-differenziali a sistemi autonomi 27

una generica equazione integro-differenziale scalare del tipo seguente:

x′(t) = f(t, x(t)) + k

∫ t

t0

1r

exp− t− s

r

g(s, x(s))ds (2.3.1)

da risolversi sull’intervallo (0,+∞). La condizione iniziale da assegnare affinche tale pro-blema sia ben posto consiste in una funzione nota x0(t) sull’intervallo (t0, 0]. Successiva-mente si introduca una variabile ausiliaria u(t) che corrisponda alla convoluzione al tempot, ovvero:

u(t) =∫ t

t0

1r

exp− t− s

r

g(s, x(s))ds (2.3.2)

Il calcolo di u′(t) fornisce:

u′(t) =g(t, x(t))

r− 1r

∫ t

t0

1r

exp− t− s

r

g(s, x(s))ds =

g(t, x(t))r

− u(t)r

(2.3.3)

Come si vede l’equazione per la u′(t) dipende solo dalla x e dalla u stessa. E pertantopossibile passare dall’equazione (2.3.1) al seguente sistema di due equazioni in due funzioniincognite:

x′(t) = f(t, x(t)) + k u(t)u′(t) = 1

r (g(t, x(t))− u(t))(2.3.4)

a cui si devono associare le condizioni iniziali x(0) = x0(0), gia necessaria per risolvere la(2.3.1), e la nuova condizione:

u(0) =∫ 0

t0

1r

expsr

g(s, x(s))ds (2.3.5)

Tale condizione e essa stessa nota, in quanto si tratta di dover calcolare analiticamente onumericamente un integrale definito di una funzione nota (infatti sull’intervallo di integra-zione x(t) = x0(t). Si sottolinea il fatto che la variabile ausiliaria introdotta coincideesattamente con il termine che tiene traccia della storia passata della soluzionedell’equazione (2.3.1). Se f e g non presentano una dipendenza diretta da t il sistema(2.3.4) risulta anche autonomo.Trattiamo ora il caso di un sistema integro differenziale di due equazioni in due incognite.Per sistemi di dimensioni maggiori il procedimento e analogo. Consideriamo per esempioil sistema: x′1(t) = f1(t, x1(t), x2(t)) + k1

∫ tt0

1r exp

− t−s

r

g1(s, x1(s), x2(s))ds

x′2(t) = f2(t, x1(t), x2(t)) + k2

∫ tt0

1r exp

− t−s

r

g2(s, x1(s), x2(s))ds

(2.3.6)

A tale sistema, da risolversi sull’intervallo (0,+∞) vanno aggiunte, come condizioni iniziali,due funzioni x10(t) e x20(t) definite sull’intervallo (t0, 0]. Procedendo in modo analogo aquanto fatto per la (2.3.1) si introducano le due funzioni ausiliarie:

u1(t) =∫ t

t0

1r

exp− t− s

r

g1(s, x1(s), x2(s))ds (2.3.7a)


u2(t) =∫ t

t0

1r

exp− t− s

r

g2(s, x1(s), x2(s))ds (2.3.7b)

le cui derivate sono:

u′1(t) =1r

(g1(t, x1(t), x2(t))−

∫ t

t0

1r

exp− t− s

r

g1(s, x1(s), x2(s))ds

)(2.3.8a)

=1r(g1(t, x1(t), x2(t))− u1(t)) (2.3.8b)

u′2(t) =1r

(g1(t, x1(t), x2(t))−

∫ t

t0

1r

exp− t− s

r

g2(s, x1(s), x2(s))ds

)(2.3.8c)

=1r(g2(t, x1(t), x2(t))− u2(t)) (2.3.8d)

In modo del tutto analogo a quanto fatto in precedenza si puo dunque passare dal sistemaintegro-differenziale (2.3.6) al sistema differenziale seguente:

x′1(t) = f1(t, x1(t), x2(t)) + k1 u1(t)u′1(t) = 1

r (g1(t, x1(t), x2(t))− u1(t))x′2(t) = f2(t, x1(t), x2(t)) + k2 u2(t)u′2(t) = 1

r (g2(t, x1(t), x2(t))− u2(t))

(2.3.9)

a cui vanno aggiunte le condizioni iniziali x1(0) = x10(0), x2(0) = x20(0) e le due seguentideterminabili a partire da x10(t) e x20(t).

u1(0) =∫ 0

t0

1r

expsr

g1(s, x1(s), x2(s))ds

u2(0) =∫ 0

t0

1r

expsr

g2(s, x1(s), x2(s))ds

Viene sottolineato ancora una volta come le funzioni ausiliarie introdotte corrispon-dano alle due storie passate delle variabili del sistema originale (2.3.6). Il proce-dimento qui illustrato verra usato nei Capitoli 3 e 4 per trasformare l’equazione integro-differenziale scalare di Verhulst-Volterra e il sistema di Lotka-Volterra con memoria insistemi differenziali autonomi.

2.4 Convergenza del nucleo alla delta di Dirac nel sensodelle distribuzioni quando r tende a 0

Come gia anticipato nella sezione 2.2, il nucleo ha l’importante ruolo di determinare lastoria passata della variabile. In particolare si e evidenziato che, qualora il nucleo mantengauna massa costante ma si concentri vicino all’origine, la storia della variabile si riduce aivalori piu recenti della stessa. A partire da tale considerazione, siamo interessati a verificarese il nucleo (2.2.1) converga, al variare di qualche parametro, alla delta di Dirac o a unsuo multiplo. Partendo dalle considerazioni fatte nella sezione 2.2 si nota agevolmente che

Convergenza del nucleo alla delta di Dirac 29

il parametro che potrebbe far tendere il nucleo in esame alla delta di Dirac o a un suomultiplo e r, il PDN. In particolare quando r viene diminuito il nucleo tende a concentrarela sua massa nell’origine pur mantenendola invariata (integrale costante e pari a k). Talecaratteristica, gia evidenziata in Figura 2.3 lascia sperare che il nucleo tenda a k voltela delta di Dirac quando il PDN tende a zero. Procedendo invece in modo rigoroso enecessario introdurre lo spazio delle funzioni test D(R). Per verificare che effettivamenteil nucleo (2.2.1) tenda a k volte la delta di Dirac nel senso delle distribuzioni occorre, inbase alla definizione A.2.1 verificare che ∀v ∈ D valga la seguente:

limr→0+

∫R

k

rexp

− t

r

H(t)v(t)dt = k

∫Rδ(t)v(t)dt = k v(0)

Si puo dimostrare la precedente procedendo come segue.

limr→0+

∫R

k

rexp

− t

r

H(t)v(t)dt = k

[lim

r→0+

∫ ∞

0

1r

exp− t

r

v(t)dt

]Integrando per parti si ottiene:

k limr→0+

[ ∫ ∞

0

1r

exp− t

r

v(t)dt

]= k lim

r→0+

[[−exp

− t

r

v(t)

]∞0

+∫ ∞

0exp

− t

r

v′(t)dt

]Si procede analizzando separatamente i contributi dei due termini contenuti fra parentesiquadre nell’ultima espressione ottenuta. Per quanto riguarda il primo di tali termini sinoti che [

− exp− t

r

v(t)

]∞0

= v(0) ∀r > 0

e di conseguenza il primo termine ha limite kv(0). Per quanto riguarda il secondo termine,considerato che v ∈ D e di conseguenza v′ limitata e nulla fuori da un compatto e cheexp

− 1

r

e anch’esso limitato, si puo applicare il Teorema A.2.1 e dedurre che tale

termine ha limite nullo. Cio dimostra che il nucleo (2.2.1) converge effettivamente nelsenso delle distribuzioni a PMN volte la delta di Dirac quando il PDN tende a zero. Sottoipotesi abbastanza blande sulla funzione con cui si convolve il nucleo vale anche il seguenterisultato:

Microteorema 2.4.1. Si consideri:

u(t) =∫ t

t0

1r

exp− t− s

r

x(s)ds. (2.4.1)

Allora se la funzione x e identificabile con una distribuzione, comunque si scelga T > 0finito si ha:

u(t) → x(t) (2.4.2)

nel senso delle distribuzioni per ogni t < T quando r → 0.

Dimostrazione. Poiche si considerano solo tempi minori di T la funzione x puo essereidentificata con la distribuzione: x(t)H(T − t)H(t − t0). Tale distribuzione ha supportocompatto e siccome u(t) puo essere vista come ((1

r exp− t−s

r

H(t))∗ (x(t)H(T − t)H(t−

t0)))(t), in base al risultato A.2.2 si conclude.

Mostriamo come tale convergenza sussista nel caso di convoluzione tra un nucleo moltoconcentrato intorno all’origine e una funzione nota. Esempi in merito sono proposti nelleFigure 2.4 e 2.5. Tuttavia tale convergenza non assicura in alcun modo che le soluzionidei modelli contenenti integrali di tipo convolutivo come quello analizzato convergano allesoluzioni dei modelli degeneri che si ottengono sostituendo all’integrale la sola funzioneconvoluta. In particolare nei capitoli dedicati all’analisi di tali modelli si tentera di verifi-care la convergenza uniforme mediante esempi numerici.Un altro caso interessante si presenta qualora la memoria sia “distribuita” su tutto ilsemiasse positivo ovvero quando r → +∞. Siamo pertanto interessati a capire comedegenerano i modelli quando r → +∞. A tal proposito vale il seguente risultato:

Microteorema 2.4.2. Per ogni funzione x non negativa e limitata e per ogni T > 0 finitosi ha che: ∫ t

t0

1r

exp− t− s

r

x(s)ds→ 0 (2.4.3)

per ogni t tale che t0 < t < T quando r → +∞.

Dimostrazione. Iniziamo col notare che se x e non negativa si ha:∫ t

t0

1r

exp− t− s

r

x(s)ds ≥ 0.

Inoltre se M = max[t0,T ](x(t)) si ha:∫ t

t0

1r

exp− t− s

r

x(s)ds ≤M

∫ t

t0

1r

exp− t− s

r

= M(1− exp

−t+ t0r

)

≤M(1− exp−T + t0

r

)

e l’ultimo membro tende a zero quando r → +∞.

Pertanto si puo pensare che comunque si scelga un tempo finale T > 0 finito il modellodegenere corrispondente sia quello senza memoria. Esempi numerici atti a verificare se visia o meno convergenza delle soluzioni dei modelli con memoria a quelli senza memoriaquando r → +∞ verranno proposti nei capitoli 3 e 4.

2.5 Esempi di convoluzioni con funzioni note

In questa parte vengono forniti alcuni esempi di convoluzione tra il nucleo (2.2.1) e al-cune funzioni note. Poiche tale convoluzione rappresenta la storia passata della variabileall’interno dell’equazione (1.3.3) e del sistema (1.6.1) si e particolarmente interessati aconvoluzioni con funzioni periodiche. Infatti la periodicita delle convoluzioni tra nucleo(2.2.1) e funzioni periodiche risulta importante qualora si pensi che esistano soluzioni pe-riodiche o cicli limite nelle (1.3.3) e/o (1.6.1). Supponiamo infatti che esista una soluzionex(t) dell’equazione di Verhulst-Volterra che tenda a un ciclo periodico. Allora anche la

Esempi di convoluzioni con funzioni note 31

Figura 2.4: Simulazione numerica della convoluzione tra la funzione sin(t)H(t) ed il nucleo (2.2.1)nel caso in cui il nucleo stesso tenda nel senso delle distribuzioni ad una delta di Dirac (r = 0.001,k = 1). Scopo di tale simulazione e la verifica numerica che la convoluzione tra il nucleo (2.2.1) ela funzione sin(t)H(t) tenda alla funzione medesima qualora il nucleo tenda alla delta di Dirac. Inordinata sono riportati i valori della convoluzione (blu) e della funzione sin(t)H(t) (arancio) mentrein ascissa il tempo.

Figura 2.5: Simulazione numerica della convoluzione tra la funzione H(t) ed il nucleo (2.2.1) nelcaso in cui il nucleo stesso tenda nel senso delle distribuzioni ad una delta di Dirac (r = 0.001,k = 1). Scopo di tale simulazione e la verifica numerica che la convoluzione tra il nucleo (2.2.1)e la funzione H(t) tenda alla funzione medesima qualora il nucleo tenda alla delta di Dirac. Inordinata sono riportati i valori della convoluzione (blu) e della funzione H(t) (arancio) mentre inascissa il tempo.


Figura 2.6: Simulazione numerica della convoluzione tra la funzione sin(t)H(t) e il nucleo (2.2.1)nel caso in cui r = k = 1. Scopo di tale simulazione e la verifica grafica del fatto che la convoluzionedel nucleo (2.2.1) con una funzione periodica tenda ad una funzione periodica dopo un periodo distabilizzazione. In ordinata sono riportati i valori della convoluzione (blu) e della funzione sin(t)H(t)(arancio) mentre in ascissa il tempo.

derivata di tale soluzione deve stabilizzarsi su un ciclo periodico. La storia della variabile sipuo dunque esprimere come differenza di termini che diventano periodici e di conseguenzadeve essa stessa tendere a un ciclo periodico. Considerazioni del tutto analoghe valgo-no relativamente al sistema di Lotka-Volterra con memoria. Non essendoci alcun risultatoteorico generale che stabilisce che la convoluzione con una funzione periodica o che diventaperiodica sia periodica o diventi tale, verifichiamo per via numerica se, nel caso di nucleo(2.2.1), cio sia possibile. Nelle Figure 2.6, 2.7, 2.8 vengono proposte le convoluzioni tra inuclei gia mostrati in Figura 2.3 e la funzione sin(ξ)H(ξ), nonche la funzione medesima.Si noti come dopo un primo periodo di stabilizzazione la convoluzione sembra assestarsisu un ciclo periodico. Questo, come peraltro gia detto, lascia intendere che almeno in certicasi sia possibile provare a cercare soluzioni o cicli periodici delle (1.3.3) e (1.6.1). Mag-giori dettagli in proposito saranno forniti nei Capitoli 3 e 4. Non si e tuttavia interessatialle sole funzioni periodiche. Un altro caso particolarmente interessante e quello in cui sieffettui la convoluzione con una funzione costante. Supponiamo infatti che la soluzionedell’equazione di Verhulst-Volterra o del sistema di Lotka-Volterra con memoria tendanoad un equilibrio. A meno di infinitesimi, da un certo tempo in poi, la soluzione dell’e-quazione sara dunque costante. Cio implica che anche la storia della variabile si debbaassestare su un valore costante per considerazioni analoghe a quelle fatte in precedenzacirca l’esistenza di cicli limite. Anche in questo caso proponiamo una Figura che mostrichiaramente come il comportamento della convoluzione sia proprio di questo tipo. La Figu-ra in questione e la 2.9. Riteniamo inoltre utile mostrare alcune convoluzioni con funzioninote, cosicche il lettore possa intuire che tipi di andamenti vengono a crearsi a seguito ditale operazione. In questo senso di seguito vengono proposte le convoluzioni tra il nucleo


Figura 2.7: Simulazione numerica della convoluzione tra la funzione sin(t)H(t) e il nucleo (2.2.1)nel caso in cui r = 2 e k = 1. Scopo di tale simulazione e la verifica grafica del fatto che laconvoluzione del nucleo (2.2.1) con una funzione periodica tenda ad una funzione periodica dopoun periodo di stabilizzazione. In ordinata sono riportati i valori della convoluzione (verde) e dellafunzione sin(t)H(t) (arancio) mentre in ascissa il tempo.

Figura 2.8: Simulazione numerica della convoluzione tra la funzione sin(t)H(t) e il nucleo (2.2.1)nel caso in cui r = 0.5 e k = 1. Scopo di tale simulazione e la verifica grafica del fatto che laconvoluzione del nucleo (2.2.1) con una funzione periodica tenda ad una funzione periodica dopoun periodo di stabilizzazione. In ordinata sono riportati i valori della convoluzione (rosso) e dellafunzione sin(t)H(t) (arancio) mentre in ascissa il tempo.


Figura 2.9: Simulazione numerica della convoluzione tra la funzione H(t) ed il nucleo (2.2.1) nelcaso in cui r = k = 1. Scopo di tale simulazione e la verifica numerica che la convoluzione trail nucleo (2.2.1) ed una costante tenda ad una costante. In ordinata sono riportati i valori dellaconvoluzione (blu) e della funzione H(t) (arancio) mentre in ascissa il tempo.

Figura 2.10: Simulazione numerica della convoluzione tra la funzione exp(t)H(t) ed il nucleo(2.2.1) nel caso in cui r = k = 1. Scopo di tale simulazione e lo studio dell’andamento dellaconvoluzione tra il nucleo (2.2.1) e una funzione nota. In ordinata sono riportati i valori dellaconvoluzione (blu) e della funzione exp(t)H(t) (arancio) mentre in ascissa il tempo.


Figura 2.11: Simulazione numerica della convoluzione tra la funzione exp(−t)H(t) ed il nucleo(2.2.1) nel caso in cui r = k = 1. Scopo di tale simulazione e lo studio dell’andamento dellaconvoluzione tra il nucleo (2.2.1) e una funzione nota. In ordinata sono riportati i valori dellaconvoluzione (blu) e della funzione exp(−t)H(t) (arancio) mentre in ascissa il tempo.

Figura 2.12: Simulazione numerica della convoluzione tra la funzione t10H(t) ed il nucleo (2.2.1)

nel caso in cui r = k = 1. Scopo di tale simulazione e lo studio dell’andamento della convoluzionetra il nucleo (2.2.1) e una funzione nota. In ordinata sono riportati i valori della convoluzione (blu)e della funzione t

10H(t) (arancio) mentre in ascissa il tempo.


prototipale (k = r = 1) e alcune funzioni note (vedasi Figure 2.10, 2.11, 2.12). Si noticome, a meno di variazioni poco significative, gli andamenti delle convoluzioni rispecchino,soprattutto a lungo termine, quelli delle funzioni convolute con il nucleo 2.2.1.

Capitolo 3

L’equazione di Verhulst-Volterra

Oggetto di questo Capitolo e l’analisi dettagliata dell’equazione integro-differenziale scalaredi Verhulst-Volterra (equazione (1.3.3)) nel caso in cui il nucleo abbia la forma (2.2.1). Inparticolare vengono proposti numerosi esempi numerici e viene fatto un confronto trai risultati ottenuti e le conclusioni del Teorema 1.3.1. Una sezione viene dedicata adanalizzare le soluzioni dell’equazione al variare del PDN, con particolare attenzione al casoin cui tenda a zero. In tale caso viene proposto un confronto con le soluzioni del modellodi controllo logistico a cui l’equazione (1.3.3) si riconduce formalmente.

3.1 Dall’equazione al sistema

Riportiamo l’equazione (1.3.3) nel caso particolare in cui K(t, s) = K(t− s) e K abbia laforma (2.2.1). In tale situazione l’equazione di Verhulst-Volterra si scrive come segue:

x′(t) = a x(t)− b x(t)2 + x(t)(k

∫ t

t0

1r

exp− t− s

r

x(s)ds

). (3.1.1)

Richiamiamo brevemente la nomenclatura:

• a > 0 prende il nome di potenziale biologico,

• b > 0 prende il nome di parametro di controllo logistico,

• r > 0 finito1 prende il nome di parametro di dispersione del nucleo (PDN),

• k > 0, k 6= b prende il nome di parametro di massa del nucleo (PMN)2.

1Il caso particolare in cui r tenda a +∞ risulta difficile da trattare analiticamente. Tale affermazionetrovera riscontro nel seguito, in particolare analizzando la matrice Jacobiana del sistema (3.1.3) associatoalla (3.1.1). Tale matrice, nel caso r = +∞ possiede un autovalore nullo e cio, dato che il sistema (3.1.3)e non lineare, impedisce di trarre delle conclusioni rigorose a partire dai metodi di linearizzazione.

2La richiesta b 6= k e legata alla particolare forma della matrice Jacobiana del sistema (3.1.3) associatoalla (3.1.1). In particolare, nel caso b = k, uno dei due equilibri del sistema non risulta ben definito. Cionon toglie che nella trattazione verra proposto qualche esempio numerico di soluzione relativa a tale casoaffinche il lettore possa farsi un’idea dei possibili andamenti.

37

38 L’equazione di Verhulst-Volterra

Si ricorda inoltre che l’equazione (3.1.1) e da risolversi sull’intervallo temporale [0,+∞).Naturalmente affinche il problema sia ben posto e necessario assegnare anche una con-dizione iniziale che, nel caso specifico, e una funzione x0(t) definita sull’intervallo [t0, 0]ove t0 e un opportuno istante temporale convenzionalmente negativo (ergo nel passato).Seguendo il ragionamento esaustivamente descritto in 2.3 e possibile introdurre una va-riabile ausiliaria e trasformare l’equazione integro-differenziale in un sistema differenzialedi due equazioni in due incognite. Scegliendo come variabile ausiliaria la storia passatadell’incognita:

u(t) =∫ t

t0

1r

exp− t− s

r

x(s)ds, (3.1.2)

si ottiene il sistema: x′(t) = x(t)

(a− b x(t) + k u(t)

)u′(t) = 1

r (x(t)− u(t))(3.1.3)

a cui vanno aggiunte le condizioni iniziali x(0) = x0(0) e:

u(0) =∫ 0

t0

1r

expsr

x0(s)ds. (3.1.4)

Come si nota facilmente tale sistema e un sistema differenziale del primo ordine non lineare,omogeneo, autonomo ed a coefficienti costanti. Esistenza e unicita locale della soluzionediscendono automaticamente dal Teorema di Cauchy-Lipschitz A.2.5. Inoltre si ricordiche nello spazio delle fasi due traiettorie non si possono mai intersecare e di conseguenzadiscende l’unicita delle soluzioni per il sistema3 in tutto l’intervallo di definizione dellestesse che puo anche arrivare a essere tutto il semiasse destro qualora la soluzione nondiverga in tempo finito. Nel seguito della trattazione il sistema verra studiato restringendoil dominio a x ≥ 0 e u ≥ 0, inoltre, salvo avviso contrario, si supporra sempre b > k.Tali ipotesi vengono introdotte per indirizzare l’analisi solo sui casi di effettivo interessemodellistico. Si noti peraltro che tali restrizioni non limitano l’interesse dell’equazione dalpunto di vista squisitamente analitico.

3.2 Analisi degli equilibri del sistema (3.1.3)

In questa sezione analizziamo gli equilibri del sistema (3.1.3) mediante le tecniche dilinearizzazione note dall’Analisi4. Il sistema presenta i seguenti punti di equilibrio:

xeq1 = ueq1 = 0; (3.2.1a)3E importante notare che l’unicita vale per il sistema. Nel seguito si analizzera principalmente solo

la variabile x(t) che e quella rilevante dal punto di vista modellistico. Per tale variabile, relativamenteall’equazione (3.1.1) non vale l’unicita in senso stretto. Infatti si possono avere, per esempio, diversesoluzioni in corrispondenza di uguale dato iniziale x(0) ma diverse funzioni x0(t). Una situazione oppostaconduce alla possibilita di avere soluzioni coincidenti con dato x(0) uguale ma x0(t) diverso (a patto chel’integrale iniziale sia pero coincidente). Inoltre potrebbe accadere che due soluzioni corrispondenti a datiiniziali diversi si intersechino. Si sottolinea comunque che fissato x(0) e x0(t) esiste una e una sola funzionex(t) soluzione dell’equazione (3.1.1) sull’intervallo di definizione. Cio e principalmente legato al fatto chel’equazione (3.1.1) e tempo dipendente. Nel seguito della trattazione, pertanto, quando si parlera di unicitasi fara sempre riferimento al sistema sebbene l’analisi sara focalizzata solo su una delle due componenti.

4Per maggiori dettagli in merito al metodo di linearizzazione si consulti l’Appendice con particolareriferimento alla sezione A.2.1.

Analisi degli equilibri del sistema (3.1.3) 39

xeq2 = ueq2 =a

b− k. (3.2.1b)

Si inizi col notare che il punto di equilibrio (3.2.1a) ha sempre interesse modellistico ecorrisponde alla particolare situazione di estinzione della specie analizzata. L’equilibrio(3.2.1b) risulta invece di interesse modellistico solo per valori di b > k. Si noti peraltro chequalora b → k+ quest’ultimo equilibrio ha coordinate che tendono all’infinito. Iniziamol’analisi della stabilita degli equilibri partendo dal primo dei due. Tale equilibrio risulta difondamentale interesse dal punto di vista modellistico. Una eventuale asintotica stabilitaporterebbe infatti all’estinzione della popolazione. Consideriamo innanzitutto la matriceJacobiana J del sistema:

J(x, u) =[a− 2bx+ ku kx

1r −1

r

]. (3.2.2)

Valutando tale matrice nelle coordinate dell’equilibrio (3.2.1a) si ottiene:

J(0, 0) =[a 01r −1

r

]. (3.2.3)

Gli autovalori associati a tale matrice sono λ1 = a e λ2 = −1r . In base al Teorema A.2.6,

ricordandosi che a, r > 0, si deduce che l’equilibrio (3.2.1a) e instabile. Trasferendo laconclusione sull’equazione (3.1.1) si ha un primo risultato importante. Tale risultato esintetizzato nel Microteorema seguente:

Microteorema 3.2.1. Si consideri l’equazione (3.1.1). Se x(0) > 0 allora la popolazionenon si estingue.

Tale risultato trova riscontro nelle simulazioni numeriche riportate nelle Figure 3.1 e 3.2.Come si puo apprezzare, anche qualora x(0) sia significativamente vicino all’equilibrio(3.2.1a), la soluzione si allontana dall’equilibrio medesimo. Valutando invece la matriceJacobiana (3.2.2) in corrispondenza dell’equilibrio (3.2.1b) si ha:

J(x, u) =[− a b

b−ka kb−k

1r −1

r

]. (3.2.4)

Essendo Tr(J) < 0 e Det(J) = ar > 0 entrambi gli autovalori della matrice sono negativi.

Cio comporta, sempre in base al Teorema A.2.6, che l’equilibrio (3.2.1b) sia localmenteasintoticamente stabile. In realta siamo riusciti a dimostrare un risultato piu forte:

Teorema 3.2.2. Si consideri l’equazione (3.1.1). Se x0(t) ≥ 0 e b > k allora la popolazionetende all’equilibrio (3.2.1b). Cio equivale a dire che il bacino di attrazione dell’equilibrio(3.2.1b) contiene (0,+∞)× (0,+∞).

Dimostrazione. Si noti innanzitutto che ogni soluzione che abbia origine nel primo qua-drante o sui semiassi positivi non puo mai uscire da [0,+∞) × [0,+∞). Infatti se x = 0anche x′ = 0 mentre u = 0 implica u′ ≥ 0 ∀x ≥ 0. Inoltre una qualunque soluzione cheabbia origine nel quadrato [0, a

b−k ]× [0, ab−k ] rimane confinata all’interno di tale quadrato.

Infatti se x = ab−k allora x′ ≤ 0 ∀ 0 ≤ u ≤ a

b−k . Fissato invece u = ab−k si ha immediata-

mente u′ ≤ 0 ∀x : 0 ≤ x ≤ ab−k . Tali considerazioni escludono l’esistenza di cicli limite.

Figura 3.1: Simulazioni numeriche relative alla soluzione dell’equazione (3.1.1). Scopo di talesimulazione e di mettere in risalto la repulsivita dell’equilibrio (3.2.1a). In particolare si e sceltacome funzione iniziale la costante x0(t) = 0.001 definita sull’intervallo [−5, 0]. In questo caso si eoptato per b = 2 > k = 1 e a = 1. Di conseguenza la soluzione del sistema rappresentata in blutende all’equilibrio (3.2.1b) rappresentato in verde. In ordinata si trova il valore della soluzionex(t) mentre in ascissa il tempo t.

Figura 3.2: Simulazioni numeriche relative alla soluzione dell’equazione (3.1.1). Scopo di talesimulazione e di mettere in risalto la repulsivita dell’equilibrio (3.2.1a). In particolare si e sceltacome funzione iniziale la costante x0(t) = 0.001 definita sull’intervallo [−5, 0]. In questo caso si eoptato per b = 0.9 < k = 1 e a = 1. Di conseguenza la soluzione del sistema rappresentata in bludiverge. In ordinata si trova il valore della soluzione x(t) mentre in ascissa il tempo t.

Figura 3.3: Simulazioni numeriche relative alla soluzione dell’equazione (3.1.1). Scopo di talesimulazione e di mettere in risalto l’attrattivita dell’equilibrio (3.2.1b). In particolare si e sceltacome funzione iniziale 2+sin(exp(−t2)) definita sull’intervallo [−5, 0], e di conseguenza x(0) > xeq.I parametri sono b = 2 > k = 1, a = 1 e r = 1. La soluzione del sistema viene rappresentata in blumentre l’equilibrio in verde. In ordinata si trova il valore della soluzione x(t) mentre in ascissa iltempo t.

Infatti e possibile dimostrare [PS2] che per sistemi bidimensionali un ciclo limite contienesempre un punto di equilibrio al suo interno. Poiche i due punti di equilibrio del sistema(3.1.3) associato alla (3.1.1) si trovano su vertici opposti del quadrato considerato e datoche una generica traiettoria non puo lasciare il quadrato una volta raggiuntolo, si deduceche non esistono cicli limite. Si vuole ora utilizzare il Teorema A.2.7 per dimostrare cheogni traiettoria che abbia un punto di origine a coordinate non negative converge all’e-quilibrio (3.2.1b) con la sola eccezione dell’origine. Consideriamo dapprima le traiettorieche originano nel quadrato [0, a

b−k ]× [0, ab−k ]. Non esistendo cicli limite all’interno di tale

regione le traiettorie devono convergere verso un punto di equilibrio e l’unico attrattivoe proprio l’equilibrio (3.2.1b). Se si ha come punto di partenza della traiettoria l’originenon si converge all’equilibrio (3.2.1b) poiche l’origine e essa stessa un punto di equilibrio.Per dimostrare il Teorema anche per generici dati iniziali di coordinate non negative siprocede notando che se u > x allora u′ < 0 e se x > u ∧ x > a

b−k si ha x′ < 0. Allora,considerata una generica traiettoria di coordinate iniziali (x0, y0) con x0 >

ab−k ∨u0 >

ab−k

tale traiettoria non potra mai uscire dal quadrato [0,max(x0, u0)]× [0,max(x0, u0)]. Nonessendoci cicli limite anche tali traiettorie convergono all’equilibrio attrattivo (3.2.1b) equesto conclude.

L’affermazione riportata nel Teorema 3.2.2 viene messa bene in risalto dalle simulazioninumeriche riportate nelle Figure 3.3 e 3.4 nonche nella gia citata Figura 3.1. Concludiamoquesta sezione trattando brevemente il caso in cui b ≤ k. Se b < k una semplice analisiteorica evidenzia che ogni soluzione rilevante dal punto di vista modellistico diverge. In-fatti l’unico equilibrio a coordinate non negative risulta repulsivo e non potendoci esserecicli limite, come gia dimostrato nel provare il Teorema 3.2.2, si conclude quanto detto.Un riscontro grafico di tale comportamento si ha nelle Figure 3.5 e 3.6. Trattiamo infine il

Figura 3.4: Simulazioni numeriche relative alla soluzione dell’equazione (3.1.1). Scopo di talesimulazione e di mettere in risalto l’attrattivita dell’equilibrio (3.2.1b). In particolare si e sceltacome funzione iniziale exp(−t2) definita sull’intervallo [−5, 0], e di conseguenza x(0) < xeq. Iparametri sono b = 2 > k = 1, a = 1.5 e r = 1. La soluzione del sistema viene rappresentata in blumentre l’equilibrio in verde. In ordinata si trova il valore della soluzione x(t) mentre in ascissa iltempo t.

caso in cui b = k. Come gia messo in rilievo in precedenza in corrispondenza a tale sceltadei parametri l’equilibrio (3.2.1b) non risulta ben definito. Di conseguenza si preferisceoptare per mostrare direttamente gli andamenti della soluzione ottenuti mediante simu-lazioni numeriche rispetto ad una analisi teorica rigorosa che risulta tecnicamente moltocomplicata. Una possibile soluzione e mostrata in Figura 3.7. A partire da tale Figurasembra che, almeno asintoticamente, la soluzione diverga in modo lineare. Non siamotuttavia riusciti a dare una dimostrazione rigorosa di tale comportamento.

Nota

In tutte le simulazioni numeriche e nelle sezioni che seguono verranno proposti solo i gra-fici relativi alla variabile x. Qualora si facciano dei commenti relativi agli andamenti dellesoluzioni dell’equazione (3.1.1), come nelle sezioni 3.3 o 3.4, si pensera sempre all’anda-mento della variabile x nel tempo ottenuto a partire dalla risoluzione numerica del sistema(3.1.3)5. Si invita dunque il lettore a prestare attenzione al fatto che nel seguito la x potraessere utilizzata sia come componente del sistema che come soluzione dell’equazione.

3.3 Confronto con il Teorema 1.3.1

Alla luce del Teorema 3.2.2, ogni soluzione corrispondente a dati iniziali di interesse mo-dellistico, qualora sia soddisfatta l’ipotesi b > k > 0 e la positivita degli altri parametri delsistema (3.1.3), tende all’equilibrio non nullo (3.2.1b). Vogliamo confrontare tale risultato

5Per risolvere numericamente il sistema si e utilizzato il solutore numerico per sistemi di ODE integratonel software Mathematica 5.2 della Wolfram Research.

Confronto con il Teorema 1.3.1 43

Figura 3.5: Simulazioni numeriche relative alla soluzione dell’equazione (3.1.1). Scopo di talesimulazione e di mettere in risalto come la soluzione diverga se b < k. In particolare si e sceltacome funzione iniziale la costante 1 definita sull’intervallo [−5, 0]. I parametri sono b = 1 < k = 1.2,a = 1 e r = 1. La soluzione del sistema viene rappresentata in blu. In ordinata si trova il valoredella soluzione x(t) mentre in ascissa il tempo t.

Figura 3.6: Simulazioni numeriche relative alla soluzione dell’equazione (3.1.1). Scopo di talesimulazione e di mettere in risalto come la soluzione diverga se b < k. In particolare si e sceltacome funzione iniziale la costante 1 definita sull’intervallo [−5, 0]. I parametri sono b = 1 < k = 2,a = 1 e r = 1. La soluzione del sistema viene rappresentata in blu. In ordinata si trova il valoredella soluzione x(t) mentre in ascissa il tempo t.

Figura 3.7: Simulazione numerica relativa alla soluzione dell’equazione (3.1.1). Scopo di talesimulazione e di mostrare il comportamento della soluzione nel particolare caso in cui b = k.L’andamento che si osserva e di tipo asintoticamente linearmente divergente. In particolare sie scelta come funzione iniziale la costante 0.5 definita sull’intervallo [−5, 0]. I parametri sonob = k = 1, a = 1 e r = 1. La soluzione del sistema viene rappresentata in blu. In ordinata si trovail valore della soluzione x(t) mentre in ascissa il tempo t.

con quanto asserito nel Teorema 1.3.1. Tale Teorema fornisce delle condizioni sufficien-ti affinche la soluzione dell’equazione (3.1.1) sia interamente contenuta all’interno di uncerto intervallo [p, P ]6 che soddisfi p < x0(t) < P ∀ t : t0 < t < 0. L’interesse di taleraffronto sara focalizzato sul valutare se, in questo caso specifico, le condizioni siano anchenecessarie e non solo sufficienti. In base al Teorema 3.2.2 e sicuramente necessario che lacoordinata x del punto di equilibrio del sistema (3.1.3) sia compresa tra p e P . Possiamopertanto impostare le due disequazioni:

a

b− k p⇒ a− bp+ kp > 0. (3.3.1b)

Considerando che k ha il significato di norma L1 del nucleo e ricordandosi le considerazionie le ipotesi precise del Teorema 1.3.1 e facile verificare che la disequazione (3.3.1a) coincidecon una delle due disequazioni facenti parte delle ipotesi del Teorema 1.3.1. Pertanto,almeno nel caso di nuclei a decadimento esponenziale del tipo 2.2.1, la condizione (3.3.1a)risulta non solo sufficiente per il Teorema 1.3.1 ma anche necessaria. Invece la disequazione(3.3.1b), per quanto necessaria al fine di ottenere un limite inferiore sul valore di equilibriodella soluzione, non trova riscontro nell’altra disequazione del Teorema 1.3.1. In particolarela disequazione facente parte delle ipotesi di tale Teorema, scritta in funzione del PMNdiventa:

a− bp− kP > 0. (3.3.2)

6Usiamo qui le lettere p e P in luogo, rispettivamente, di r e R usate nel Teorema 1.3.1 per evitarepossibili ambiguita tra tali soglie e il PDN.

Figura 3.8: Simulazioni numeriche relative alla soluzione dell’equazione (3.1.1). Lo scopo di talesimulazione e quello di mostrare che il fatto che l’equilibrio (3.2.1b) sia maggiore di p e condizionesolo necessaria ma non sufficiente affinche x(t) > p ∀t > 0. Infatti, scelto p = xeq, non e comunqueverificato x(t) > p ∀t > 0 anche se la soluzione tende all’equilibrio (3.2.1b). I parametri con cui sie effettuata tale simulazione sono a = 1, b = 2, k = 1, r = 10 mentre la funzione iniziale e x0(t) = 2definita su [−5, 0]. In blu viene riportata la soluzione x(t) mentre in verde l’equilibrio (3.2.1b).

Le ragioni di tale mancata corrispondenza sono duplici. In primo luogo bisogna pensareche la disequazione (3.3.2) e solo sufficiente e non necessaria e pertanto piu restrittiva.Inoltre la condizione dettata dalla (3.3.1b) e necessaria affinche il valore di equilibrio siamaggiore della soglia p ma non fornisce alcuna garanzia sul fatto che una soluzione cheparta sopra la soglia non possa scendere, avvicinandosi all’equilibrio, al di sotto di talevalore. In effetti, come evidenziato per esempio nelle Figure 3.8, 3.9, 3.10 esistono delleparticolari combinazioni dei parametri che forniscono gli andamenti appena citati. Etuttavia possibile notare che sembra esserci un limite inferiore oltre il quale soluzioniche partono da un valore x0(0) = x(0) > a

b−k non possano scendere. Dalle simulazioninumeriche sembra che tale valore limite sia pari ad a

b . In realta tale risultato puo esseredimostrato e viene racchiuso nel seguente Microteorema:

Microteorema 3.3.1. Si consideri l’equazione (3.1.1). Se b > k e x0(0) > ab−k allora la

popolazione non puo mai scendere sotto la soglia di ab individui.

Dimostrazione. Supponiamo per assurdo che la x attraversi la retta di equazione x = ab .

Allora deve essere x′ < 0. Ma se x = ab si ha x′ > 0: assurdo.

Si noti che la soglia oltre la quale la popolazione non scende mai, qualora parta da un datoiniziale x(0) > a

b−k , coincide esattamente col valore dell’equilibrio nel caso in cui non vi siala memoria. Ovvero se il PDN tende a infinito (nucleo disperso su tutto il semiasse positivoma con valori infinitesimi) l’equazione sembra tendere all’equazione senza memoria. Talerisultato puo essere formalizzato nel caso b > k.

Teorema 3.3.2. Se b > k la soluzione del modello (3.1.1) e sempre maggiore o uguale diquella del modello di controllo logistico (1.2.1) per ogni t tale che t0 ≤ t ≤ T ove T e un


Figura 3.9: Simulazioni numeriche relative alla soluzione dell’equazione (3.1.1). Lo scopo di talesimulazione e quello di mostrare che il fatto che l’equilibrio (3.2.1b) sia maggiore di p e condizionesolo necessaria ma non sufficiente affinche x(t) > p ∀t > 0. Infatti, scelto p = xeq, non e comunqueverificato x(t) > p ∀t > 0 anche se la soluzione tende all’equilibrio (3.2.1b). Inoltre si vuolemostrare che esiste un valore minimo sotto il quale la soluzione non puo scendere se x(0) > xeq.Tale valore e a

b . I parametri con cui si e effettuata tale simulazione sono a = 1, b = 2, k = 1, r = 50mentre la funzione iniziale e x0(t) = 2 definita su [−5, 0]. In blu viene riportata la soluzione x(t)mentre in verde l’equilibrio (3.2.1b).

Figura 3.10: Simulazioni numeriche relative alla soluzione dell’equazione (3.1.1). Lo scopo di talesimulazione e quello di mostrare che il fatto che l’equilibrio (3.2.1b) sia maggiore di p e condizionesolo necessaria ma non sufficiente affinche x(t) > p ∀t > 0. Infatti, scelto p = xeq, non e comunqueverificato x(t) > p ∀t > 0 anche se la soluzione tende all’equilibrio (3.2.1b). Inoltre si vuolemostrare che esiste un valore minimo sotto il quale la soluzione non puo scendere se x(0) > xeq.Tale valore e a

b . I parametri con cui si e effettuata tale simulazione sono a = 1, b = 2, k = 1,r = 250 mentre la funzione iniziale e x0(t) = 2 definita su [−5, 0]. In blu viene riportata la soluzionex(t) mentre in verde l’equilibrio (3.2.1b).

opportuno istante temporale positivo finito. Inoltre la soluzione della (3.1.1) tende a quelladel corrispondente modello di controllo logistico per ogni t tale che t0 < t < T quando ilPDN tende a infinito.

Dimostrazione. Consideriamo le soluzioni x dell’equazione di Verhulst-Volterra e y delcorrispondente modello senza memoria y′ = ay − by2 con uguale dato iniziale x(t0) =y(t0). Ricordiamo che in base alla soluzione analitica (1.2.2) del modello senza memoriae ai risultati proposti in questo capitolo x e y sono limitate e non negative. Sia dunquec = min[t0,T ](x(t), y(t)) e C = max[t0,T ](x(t), y(t)). Essendo x e y soluzioni della (3.1.1) edel relativo modello senza memoria dovra essere:

x′ − y′ = a(x− y)− b(x− y)(x+ y) + kx

∫ t

t0

1r

exp− t− s

r

x(s)ds.

Chiamando w la differenza tra x e y si ha:

w′ = aw − bw(x+ y) + kx

∫ t

t0

1r

exp− t− s

r

x(s)ds

con w(t0) = 0. Poiche x e sovrasoluzione dell’equazione y′ = ay − by2, grazie ai criteri diconfronto per soluzioni di EDO si ha w ≥ 0. Notiamo ora che per t0 < t < T si ha:

0 ≤ kx

∫ t

t0

1r

exp− t− s

r

x(s)ds ≤ kC2(1− exp

− T − t0

r

) = ε(r)7.

Possiamo ora maggiorare e minorare w′ nel seguente modo:

(a− bC)w ≤ w′ ≤ (a− bc)w + ε(r).

Risolvendo le due equazioni differenziali relative ai limiti della precedente disuguaglianza,sempre grazie ai criteri di confronto per soluzioni di EDO si ha che per t0 < t < T :

0 ≤ w(t) ≤(−1 + exp

(a− bc)(t− t0)

)

a− bcε(r).

Tenendo conto che l’ultimo termine, qualora sia ben definito8, e sempre positivo e tendentea zero quando r tende a infinito si conclude. In particolare fissato un qualunque T > 0finito si ha che x(t) ≥ y(t) per ogni t0 < t < T . Inoltre comunque scelta una costanteδ > 0 esiste r tale che x− y ≤ δ per ogni t tale che t0 < t < T .

Esempi significativi che mostrano la validita del Teorema 3.3.2 sono mostrati nelle Figure3.11 e 3.12.Alla luce dei risultati esposti si nota facilmente che una condizione necessaria affinche sia

7Ove con ε(r). si e indicata una quantita che tende a zero quando r → +∞. Infatti kC2(1 − expn−

T−t0r

o) tende a zero quando r tende a infinito.

8L’unico caso in cui tale termine e non definito e quello in cui a − bc = 0 ovvero c = ab. Qualora

tale situazione si verifichi si puo dimostrare la convergenza scegliendo nella dimostrazione in luogo di c unarbitrario valore c positivo e minore di a

b.


Figura 3.11: Simulazioni numeriche della soluzione dell’equazione (3.1.1) blu e del corrispettivomodello senza memoria magenta. Scopo di tale simulazione e quello di mostrare come la soluzionedella (3.1.1) tenda a quella del modello degenere quando r tende a infinito. Per i parametri si sonoscelti i seguenti valori: a = 1, b = 1.5, k = 1, r = 2000. Come dato iniziale si e assegnata lacostante x0 = 2. In ordinata e riportato il valore numerico delle due soluzione mentre in ascissa iltempo.

Figura 3.12: Simulazioni numeriche della soluzione dell’equazione (3.1.1) blu e del corrispettivomodello senza memoria magenta. Scopo di tale simulazione e quello di mostrare come la soluzionedella (3.1.1) tenda a quella del modello degenere quando r tende a infinito. Per i parametri si sonoscelti i seguenti valori: a = 1, b = 1.5, k = 1, r = 5000. Come dato iniziale si e assegnata lacostante x0 = 0.5. In ordinata e riportato il valore numerico delle due soluzione mentre in ascissail tempo.

Soluzioni al variare del PDN 49

valido l’asserto del Teorema 1.3.1 e la seguente:

a− bp > 0. (3.3.3)

Tale condizione e meno restrittiva della (3.3.2) come del resto ci si poteva aspettare.La (3.3.2) risulta infatti solo sufficiente mentre la (3.3.3) e una condizione necessaria esufficiente. Sintetizziamo i risultati ottenuti nel seguente Teorema la cui dimostrazioneviene omessa in quanto banale conseguenza dei risultati gia proposti.

Teorema 3.3.3. Si consideri l’equazione (3.1.1) nel caso in cui k < b. Siano p, P costantipositive tali che per t0 ≤ s ≤ 0 si abbia 0 0.

Esempio 3.3.1. Verifichiamo la validita di tale Teorema nel caso in cui la (3.1.1) degenerinel seguente modello di controllo logistico senza memoria (per esempio qualora k = 0 or = +∞):

x′(t) = a x(t)− b x2(t)

La soluzione di tale equazione e stata trattata in dettaglio nella sezione 1.2. In particolarese x0 < xeq = a

b la soluzione e monotona crescente e tende all’equilibrio mentre se x(0) >xeq la soluzione e monotona decrescente e tende all’equilibrio. Le due disuguaglianze delTeorema 3.3.3 si riducono in questo caso degenere a:

a− bP < 0 ⇒ a

b 0 ⇒ a

b> p

ove x(0) ∈ (p, P ). Il Teorema afferma pertanto che p < x(t) 0 se e solo sel’equilibrio e il dato iniziale sono compresi tra p e P e questo trova riscontro con quantodetto sulle soluzioni del modello senza memoria.

3.4 Soluzioni al variare del PDN

In questa ultima sezione si analizzano le caratteristiche salienti che vengono evidenziatedalle simulazioni numeriche al variare del parametro r e, ove possibile, si fornisce una giu-stificazione teorica che corrobori l’aspetto numerico. Iniziamo col sottolineare l’andamentodegli autovalori delle matrici Jacobiane valutate negli equilibri del sistema al variare delparametro r. In particolare siamo interessati al caso in cui b > k. Nelle figure 3.13 e 3.14sono illustrati gli andamenti degli autovalori di tali matrici. Come si nota, l’autovalore po-sitivo della matrice Jacobiana associata all’equilibrio nullo (3.2.1a) rimane costante mentrequello negativo tende a zero al crescere di r. Un discorso analogo vale per l’equilibrio nonnullo (3.2.1b). In questo caso tuttavia ambedue gli autovalori sono negativi e crescenti al


Figura 3.13: Simulazioni numeriche relative agli andamenti degli autovalori della matrice Jacobia-na associata all’equilibrio (3.2.1a). Lo scopo di tale simulazione e quello di mostrare il fatto che unodi tali autovalori, rappresentato in rosso rimane costante e pari ad a mentre l’altro, rappresentatoin blu tende a zero. Per i parametri si sono scelti i seguenti valori: a = 1, b = 2, k = 1. In ordinatae riportato il valore numerico degli autovalori mentre in ascissa il valore di r corrispondente.

crescere di r: uno tende a zero mentre l’altro tende ad un valore costante. Si vuole enfatiz-zare il fatto che c’e evidenza numerica per concludere che tali andamenti, nel caso b > k,non dipendono dalla particolare scelta dei parametri se non per il valore della costante nonnulla cui tende un autovalore associato all’equilibrio (3.2.1b) e per quello dell’autovalorecostante associato all’equilibrio nullo (3.2.1a). Essendo il modulo degli autovalori dei dueequilibri legato alla velocita di attrazione o repulsione dei medesimi si puo concludere che:

al crescere di r diminuisce la velocita con cui le soluzioni tendono all’equilibrio nonnullo (3.2.1b) e con cui soluzioni che partono vicino all’equilibrio nullo (3.2.1a) siallontanano dal medesimo.

al diminuire di r aumenta la velocita con cui le soluzioni tendono all’equilibrio (3.2.1b)e con cui soluzioni che partono vicino all’equilibrio (3.2.1a) si allontanano dal mede-simo.

Tali comportamenti sono chiaramente mostrate nella Figura 3.15. Ci soffermiamo orabrevemente sull’analisi del caso particolarmente significativo in cui r → 0+. Come giadiscusso nella sezione 2.4 in tale caso il nucleo (2.2.1) tende, nel senso delle distribuzioni, ak volte la distribuzione δ di Dirac. Ricordandosi che tale distribuzione e l’elemento neutrodel prodotto di convoluzione, sostituendo formalmente la stessa all’interno dell’equazione(3.1.1) si ottiene:

x′(t) = x(t)(a− (b− k)x(t)). (3.4.1)

La (3.4.1) rappresenta, per valori di b > k, l’equazione del modello biologico di controllologistico in cui il parametro a ha ruolo di potenziale biologico mentre il termine a

b−krappresenta la capacita dell’ambiente. Siamo pertanto interessati a vedere se la soluzionedell’equazione (3.1.1) tenda alla soluzione dell’equazione (3.4.1). Una verifica rigorosa


Figura 3.14: Simulazioni numeriche relative agli andamenti degli autovalori della matrice Jaco-biana associata all’equilibrio (3.2.1a). Lo scopo di tale simulazione e quello di mostrare il fatto cheuno di tali autovalori, rappresentato in rosso tende a zero mentre l’altro, rappresentato in blu tendead una costante. Per i parametri si sono scelti i seguenti valori: a = 1, b = 2, k = 1. In ordinata eriportato il valore numerico degli autovalori mentre in ascissa il valore di r corrispondente.

Figura 3.15: Simulazioni numeriche relative alla soluzione dell’equazione (3.1.1). Scopo di taleFigura e quello di mettere in risalto le diverse velocita di repulsione/attrazione dei due equilibri(3.2.1a) ed (3.2.1b) al variare del PDN. In particolare maggiore e il parametro r minore e la velocitadi attrazione/repulsione dei due equilibri. In verde la soluzione per r = 0.1, in rosso per r = 1 edin blu per r = 10. In nero sono riportati l’equilibrio (3.2.1b) (tratteggiato) e il dato iniziale x0(t)definito sull’intervallo [−5, 0]. Gli altri parametri relativi alla simulazione sono a = 1, b = 2, k = 1.In ordinata e riportato il numero di individui mentre in ascissa il tempo.


Figura 3.16: Simulazioni numeriche della soluzione dell’equazione (3.1.1) blu e (3.4.1) magenta.Scopo di tale simulazione e quello di mostrare come la soluzione della (3.1.1) tenda a quella della(3.4.1) quando r tende a zero. Per i parametri si sono scelti i seguenti valori: a = 1, b = 1.5, k = 1,r = 0.02. Come dato iniziale si e assegnata la costante x0 = 1. In ordinata e riportato il valorenumerico delle due soluzione mentre in ascissa il tempo.

Figura 3.17: Simulazioni numeriche della soluzione dell’equazione (3.1.1) blu e (3.4.1) magenta.Scopo di tale simulazione e quello di mostrare come la soluzione della (3.1.1) tenda a quella della3.4.1 quando r tende a zero. Per i parametri si sono scelti i seguenti valori: a = 1, b = 1.5, k = 1,r = 0.02. Come dato iniziale si e assegnata la costante x0 = 3. In ordinata e riportato il valorenumerico delle due soluzione mentre in ascissa il tempo.


richiede l’introduzione di opportuni spazi funzionali e la scelta di una opportuna norma.Stime di convergenza si possono effettuare mediante tecniche che esulano dagli scopi diquesta trattazione. Ci accontentiamo pertanto di mostrare come dal punto di vista graficosembra esserci convergenza uniforme9 tra le soluzioni dei modelli. Si noti, per inciso, chela convergenza uniforme e una convergenza molto forte che implica le convergenze nellenorme integrali. In particolare si faccia riferimento alle Figure 3.16 e 3.17, notando cheil parametro r e stato scelto sufficientemente grande affinche si potesse avere una minimadistinzione tra i grafici delle soluzioni della (3.1.1) e (3.4.1). Diminuendo ulteriormente ilPDN i due grafici risultano infatti indistinguibili.

9Tale affermazione e riferita al caso b > k. Nel caso b = k, anche per valori di r molto piccoli, sembrache il comportamento della soluzione dell’equazione di Verhulst-Volterra sia di tipo divergente in modolineare discostandosi significativamente dal corrispondente modello degenere che prevede una soluzionecon crescita esponenziale.

Capitolo 4

Il sistema di Lotka-Volterra conmemoria

Oggetto di questo Capitolo e l’analisi dettagliata del sistema integro-differenziale di Lotka-Volterra con memoria (equazione (1.6.1)) nel caso in cui il nucleo abbia la solita formaespressa nella (2.2.1). In particolare vengono proposti numerosi esempi numerici e qualoraquesti ci abbiano spinto a formulare delle ipotesi, tali ipotesi sono formalizzate in Teoremio Congetture. Una sezione viene dedicata ad analizzare le soluzioni del sistema al variaredei PDN, con particolare attenzione al caso in cui tendano a zero o a infinito. In tali casivengono proposti i confronti con le soluzioni dei modelli di Lotka-Volterra con controllologistico a cui il sistema (1.6.1) si riconduce formalmente. Inoltre si analizza l’esistenza dieventuali orbite periodiche e/o cicli limite.

4.1 Dal sistema integro-differenziale a quello autonomo

Riportiamo il sistema (1.6.1) nel caso particolare in cui K(t, s) = K(t − s) e K abbia laforma (2.2.1). In tale situazione il sistema di Lotka-Volterra con memoria si scrive: x′(t) = x(t)

(a− bx(t)− gy(t)− k1

∫ tt0

1r1

exp− t−s

r1

y(s)ds

)y′(t) = y(t)

(− α− βy(t) + γx(t) + k2

∫ tt0

1r2

exp− t−s

r2

x(s)ds

) . (4.1.1)

Richiamiamo brevemente la nomenclatura:

• x(t) e la variabile indicante il numero di prede, y(t) di predatori,

• a > 0 e α > 0 prendono il nome di potenziali biologici,

• b ≥ 0 e β ≥ 0 prendono il nome di parametri di controllo logistico,

• g ≥ 0 e γ ≥ 0 prendono il nome di parametri di interazione istantanea preda-predatore,

• k1 ≥ 0, k2 ≥ 0 prendono il nome di parametri di massa dei nuclei (PMN),

• r1 > 0, r2 > 0 finiti prendono il nome di parametri di dispersione dei nuclei (PDN).

55

56 Il sistema di Lotka-Volterra con memoria

Si richiama inoltre che il sistema (4.1.1) e da risolversi sull’intervallo temporale [0,+∞).Affinche il problema sia ben posto e necessario assegnare due funzioni iniziali x0(t) ey0(t) sull’intervallo [t0, 0], ove t0 e un opportuno istante temporale convenzionalmentenegativo. Seguendo il procedimento esaustivamente descritto in 2.3 e possibile introdurredue variabili ausiliarie u e v e trasformare il sistema integro-differenziale in un sistemadifferenziale autonomo di quattro equazioni in quattro funzioni incognite. Scegliendo comevariabili ausiliarie le seguenti:

u(t) =∫ t

t0

1r2

exp− t− s

r2

x(s)ds; (4.1.2a)

v(t) =∫ t

t0

1r1

exp− t− s

r1

y(s)ds. (4.1.2b)

si perviene al sistema:

x′(t) = x(t)(a− bx(t)− gy(t)− k1v(t)

)y′(t) = y(t)

(− α− βy(t) + γx(t) + k2u(t)

)v′(t) = 1

r1

(y(t)− v(t)

)u′(t) = 1

r2

(x(t)− u(t)

) , (4.1.3)

a cui vanno aggiunte le condizioni iniziali x(0) = x0(0), y(0) = y0(0) e le due seguentideterminabili a partire dalle funzioni x0(t) e y0(t):

u(0) =∫ 0

t0

1r2

exp sr2

x0(s)ds,

v(0) =∫ 0

t0

1r1

exp sr1

y0(s)ds.

Viene sottolineato che le variabili ausiliarie corrispondono alle storie passate delle variabilioriginali del sistema. Infine notiamo ancora una volta come il sistema (4.1.3) sia non linearema autonomo e pertanto i suoi equilibri si possono analizzare col metodo di linearizzazionedescritto in appendice (A.2.1). Esistenza e unicita locale delle soluzioni del sistema (4.1.3)sono garantite dal Teorema di Cauchy-Lipschitz A.2.5. Inoltre si ricordi che nello spaziodelle fasi di un sistema autonomo due traiettorie non si possono mai intersecare e di con-seguenza si ha l’unicita1delle componenti della soluzione su tutto l’intervallo di definizioneche puo coincidere con tutto il semiasse destro qualora le medesime non esplodano prima.Nel seguito della trattazione il sistema viene studiato restringendo il dominio a x ≥ 0,y ≥ 0, v ≥ 0 e u ≥ 0. Tale ipotesi viene introdotta per indirizzare l’analisi solo sui casi diinteresse modellistico.

1E importante notare che l’unicita vale per il sistema. Nel seguito si fara riferimento essenzialmentealle sole componenti x(t) e y(t) che sono quelle rilevanti dal punto di vista modellistico. Per tali variabili,relativamente al sistema (4.1.1) non vale l’unicita in senso stretto. Infatti si possono avere diverse soluzionix(t) o y(t) in corrispondenza di uguali dati iniziali x(0) e/o y(0) ma diverse funzioni iniziali x0(t) e y0(t). Einoltre possibile che vi siano soluzioni x(t) e/o y(t) coincidenti ma corrispondenti a diverse funzioni inizialix0(t) e y0(t) (a patto che gli integrali iniziali, x(0) e y(0) coincidano). Queste sono solo alcune dei possibilicasi che si possono presentare e che devono mettere in guardia il lettore dal confondere l’unicita per lesoluzioni del sistema (4.1.3) e quelle del sistema (4.1.1).


4.2 Analisi degli equilibri del sistema (4.1.3)

In questa sezione analizziamo gli equilibri del sistema (4.1.3) mediante le tecniche dilinearizzazione note dall’Analisi. Il sistema presenta i seguenti punti di equilibrio:

xeq1 = ueq1 = 0,yeq1 = veq1 = 0; (4.2.1a)

xeq2 = ueq2 =a

b,

yeq2 = veq2 = 0; (4.2.1b)

xeq3 = ueq3 = 0,

yeq3 = veq3 = −αβ

; (4.2.1c)

xeq4 = ueq4 =aβ + αk1 + αg

gk2 + γk1 + gγ + bβ + k1k2,

yeq4 = veq4 =aγ + ak2 − bα

gk2 + γk1 + gγ + bβ + k1k2. (4.2.1d)

Per prima cosa notiamo che l’equilibrio (4.2.1c) ha sempre una coordinata negativa seβ > 0 (o non e definito qualora β = 0). Di conseguenza, tale equilibrio non ha alcunarilevanza dal punto di vista modellistico anche considerato che soluzioni con dati iniziali acoordinate non negative permangono non negative per ogni tempo t > 0. Tale affermazionee banale conseguenza del fatto che, considerando il sistema (4.1.1), se x = 0 anche x′ = 0e se y = 0 anche y′ = 0. L’analisi viene pertanto limitata ai soli equilibri (4.2.1a), (4.2.1b)ed (4.2.1d). Notiamo che tutti gli equilibri coincidono con quelli del modellosenza memoria (1.5.1) se si pone k1 = k2 = 0. Questa osservazione e importantein quanto mostra come il sistema considerato sia un’estensione di quello senzamemoria e che, qualora i termini di memoria siano evanescenti, gli equilibridei due modelli coincidono. Prima di passare a tale analisi e tuttavia convenienteintrodurre la matrice Jacobiana associata al sistema (4.1.3) nella sua forma piu generale(ovvero come funzione delle coordinate):

J(x, y, v, u) =

a− k1v − 2bx− gy −gx −k1x 0

γy k2u− α− 2βy + γx 0 k2y0 1

r1− 1

r10

1r2

0 0 − 1r2

. (4.2.2)

4.2.1 Analisi dell’equilibrio nullo (4.2.1a)

L’equilibrio (4.2.1a) e molto significativo dal punto di vista modellistico. Se fosse local-mente asintoticamente stabile infatti esisterebbero delle particolari coppie di dati iniziali


Figura 4.1: Soluzione numerica del sistema (4.1.1) atta a mostrare la repulsivita dell’equilibrio(4.2.1a). Si noti come partendo anche molto vicino all’equilibrio in questione (x0(t) = y0(t) =0.001) le componenti della soluzione si allontanano dal valore nullo. In particolare in ordinatasono mostrati gli andamenti del numero di prede (rosso) e di predatori (blu) mentre in ascissa iltempo. I parametri relativi a questa simulazione sono: a = 0.5, b = 1, g = 0, α = 1, β = 0.5,γ = 2, k1 = k2 = r1 = r2 = 1. Infine poiche in corrispondenza a tale scelta dei parametril’equilibrio (4.2.1b) (di coordinate xeq2 = 0.5 e yeq2 = 0) risulta repulsivo la soluzione non tende atale equilibrio.

per le quali entrambe le popolazioni si estinguerebbero. Per stabilire la natura di talepunto di equilibrio usiamo il metodo di linearizzazione illustrato in A.2.1. Valutando lamatrice Jacobiana (4.2.2) nell’equilibrio (4.2.1a) si ha:

J(0, 0, 0, 0) =

a 0 0 00 −α 0 00 1

r1− 1

r10

1r2

0 0 − 1r2

. (4.2.3)

Gli autovalori di tale matrice, calcolati mediante il Software simbolico Mathematica sono:λ1 = a, λ2 = − 1

r1, λ3 = − 1

r2e λ4 = −α. Come si nota, in base alle ipotesi fatte sui

parametri, esiste sempre un autovalore con parte reale positiva (λ1). Pertanto, applicandoil Teorema A.2.6 si conclude che l’equilibrio (4.2.1a) e instabile. Tale affermazione trovariscontro nelle Figure 4.1 (x e y in funzione del tempo) e 4.2 (traiettorie nel piano xy)che illustrano come una soluzione relativa a dati iniziali particolarmente vicini all’originesi discosti dall’origine stessa.

4.2.2 Analisi dell’equilibrio (4.2.1b) (estinzione dei predatori)

Dal punto di vista modellistico, l’equilibrio (4.2.1b) corrisponde alla particolare situazionein cui una delle due popolazioni (predatori) si estingue mentre l’altra (prede) si stabilizzasu una quantita prefissata e pari a a

b . Si noti che tale equilibrio e ben definito solo se epresente il termine di controllo logistico sulla specie predata. Inoltre e interessante notarecome tale equilibrio coincida con l’equilibrio (1.5.2b) del sistema di Lotka-Volterra con

Analisi dell’equilibrio (4.2.1b) (estinzione dei predatori) 59

Figura 4.2: Soluzione numerica del sistema (4.1.1) atta a mostrare la repulsivita dell’equilibrio(4.2.1a). Si noti come partendo anche molto vicino all’equilibrio in questione (x0(t) = y0(t) =0.001) le componenti della soluzione si allontanano dal valore nullo. In particolare viene mostratala traiettoria della soluzione (verde) nel piano xy. In ordinata e presente il numero di predatorimentre in ascissa il numero di prede. I parametri relativi a questa simulazione sono: a = 0.5, b = 1,g = 0, α = 1, β = 0.5, γ = 2, k1 = k2 = r1 = r2 = 1. Infine poiche in corrispondenza a talescelta dei parametri l’equilibrio (4.2.1b) (di coordinate xeq2 = 0.5 e yeq2 = 0) risulta repulsivo lasoluzione non tende a tale equilibrio.

controllo logistico ma senza memoria. Per condurre tale analisi si utilizza la solita tecnicadi linearizzazione. La matrice Jacobiana (4.2.2) del sistema valutata in corrispondenzadell’equilibrio (4.2.1b) assume la seguente forma:

J(a

b, 0, 0,

a

b) =

−a −ag

b −ak1b 0

0 k2a+aγ−bαb 0 0

0 1r1

− 1r1

01r2

0 0 − 1r2

. (4.2.4)

Gli autovalori associati a tale matrice sono: λ1 = −a, λ2 = − 1r1

, λ3 = − 1r2

e λ4 =k2a+aγ−bα

b . Come si nota tre autovalori hanno sempre parte reale negativa mentre il quartoha parte reale negativa se e solo se bα > k2a+aγ. Si osservi che tale richiesta equivale a direche l’equilibrio (4.2.1b) ha tutti gli autovalori negativi se e solo se l’equilibrio (4.2.1d) none interessante dal punto di vista modellistico in quanto possiede una coordinata negativa.Cio implica il seguente:

Microteorema 4.2.1. Se l’equilibrio (4.2.1d) (coesistenza di prede e predatori) non harilevanza dal punto di vista modellistico l’equilibrio (4.2.1b) (estinzione dei predatori) elocalmente asintoticamente stabile, viceversa e instabile.

Sarebbe interessante provare a verificare se l’equilibrio (4.2.1b), nel caso in cui sia attrat-tivo, sia non solo localmente ma anche globalmente asintoticamente stabile. Un risultatoimportante che garantisce tale proprieta viene enunciato piu avanti nel Teorema 4.3.3.Proviamo a interpretare la condizione sotto la quale l’equilibrio estinzione del predatore

Figura 4.3: Soluzione numerica del sistema (4.1.1) atta a mostrare l’attrattivita dell’equilibrio(4.2.1b). Come funzioni iniziali si sono scelte le costanti x0(t) = y0(t) = 0.5. In particolare inordinata sono mostrati gli andamenti del numero di prede (rosso) e di predatori (blu) mentre inascissa il tempo. I parametri relativi a questa simulazione sono: a = 1, b = 1, g = 0, α = 1,β = 0.5, γ = 0.5, k1 = 1, k2 = 0.5, r1 = r2 = 1.

risulta attrattivo in modo da metterne in risalto il significato modellistico. Riscriviamotale condizione nella forma:

a

b<

α

k2 + γ.

Notiamo che ab e il valore a cui tende la popolazione di prede se l’equilibrio e attratti-

vo ovvero la capacita dell’ambiente per le prede. Pertanto la condizione significa che ipredatori si estinguono se la capacita dell’ambiente relativa alle prede e suf-ficientemente piccola. Tale condizione e del resto piuttosto intuitiva: se non ci sonoabbastanza prede i predatori non riescono a sopravvivere. Simulazioni numeriche relativea soluzioni che tendono all’equilibrio (4.2.1b) sono proposte nelle Figure 4.3, 4.4, 4.5, 4.6.Inoltre nelle Figure 4.1 e 4.2 si nota come qualora l’equilibrio (4.2.1d) abbia entrambele coordinate positive la soluzione non tende all’equilibrio (4.2.1b) che per le due Figurecitate ha coordinate: xeq2 = ueq2 = 0.5, yeq2 = veq2 = 0.

4.2.3 Analisi dell’equilibrio (4.2.1d)(coesistenza di entrambe le specie)

Dal punto di vista modellistico, l’equilibrio (4.2.1d) corrisponde alla particolare situazionein cui entrambe le popolazioni si stabilizzano su una quantita prefissata. Tale equilibrio,sotto le ipotesi fatte, risulta ben definito nella maggior parte dei casi.2 Si noti inoltre chetale equilibrio ha interesse modellistico solo se bα e minore di k2a + aγ. Invitiamo inol-tre il lettore a notare che qualora si verifichi bα = k2a + aγ l’equilibrio (4.2.1d) coincideesattamente con l’equilibrio (4.2.1b). Si e interessati ad analizzare la stabilita dell’equi-librio (4.2.1d). L’idea e di sfruttare la consueta tecnica di linearizzazione. Purtroppo lamatrice Jacobiana (4.2.2) del sistema valutata in corrispondenza dell’equilibrio (4.2.1d)non si semplifica in modo significativo e anche mediante i software di calcolo simbolico

2Basta infatti che almeno uno dei seguenti prodotti sia non nullo: gk2, γk1, gγ, bβ, k1k2.

Analisi dell’equilibrio (4.2.1d)(coesistenza di entrambe le specie) 61

Figura 4.4: Soluzione numerica del sistema (4.1.1) atta a mostrare l’attrattivita dell’equilibrio(4.2.1b). Come funzioni iniziali si sono scelte le costanti x0(t) = y0(t) = 0.5. In particolare vienemostrata la traiettoria della soluzione nel piano xy. In ascissa e presente il numero di prede mentrein ordinata il numero di predatori. I parametri relativi a questa simulazione sono: a = 1, b = 1,g = 0, α = 1, β = 0.5, γ = 0.5, k1 = 1, k2 = 0.5, r1 = r2 = 1.

Figura 4.5: Soluzione numerica del sistema (4.1.1) atta a mostrare l’attrattivita dell’equilibrio(4.2.1b). Come funzioni iniziali si sono scelte le costanti x0(t) = 5 e y0(t) = 3. In particolare inordinata sono mostrati gli andamenti del numero di prede (rosso) e di predatori (blu) mentre inascissa il tempo. I parametri relativi a questa simulazione sono: a = 1, b = 1, g = 0, α = 1,β = 0.5, γ = 0.5, k1 = 1, k2 = 0.5, r1 = r2 = 1.


Figura 4.6: Soluzione numerica del sistema (4.1.1) atta a mostrare l’attrattivita dell’equilibrio(4.2.1b). Come funzioni iniziali si sono scelte le costanti x0(t) = 5 e y0(t) = 3. In particolare vienemostrata la traiettoria della soluzione nel piano xy. In ascissa e presente il numero di prede mentrein ordinata il numero di predatori. I parametri relativi a questa simulazione sono: a = 1, b = 1,g = 0, α = 1, β = 0.5, γ = 0.5, k1 = 1, k2 = 0.5, r1 = r2 = 1.

Mathematica 5.2 o Maple 10 non si riescono ad ottenere gli autovalori del sistema in for-ma chiusa. Pertanto evitiamo di riportare tale matrice. Per non lasciare la trattazionein sospeso sulla natura dell’equilibrio analizzato mostriamo i grafici degli autovalori dellamatrice Jacobiana (4.2.2) associata al sistema (4.1.3) al variare dei PDN e mantenendofissi gli altri parametri in due casi particolari:

• Nel primo caso, quello relativo a Figura 4.7 tutti gli autovalori sono a parte realenegativa per ogni valore di r1 e r2. Per quanto i grafici riportati siano relativi ad unsolo esempio, le numerose prove numeriche effettuate consentono di affermare che incerti casi l’equilibrio (4.2.1d) e attrattivo indipendentemente dalla scelta di r1 e r2.Tale situazione e analoga a quella trattata nell’esempio di sezione 4.2.4 in assenzadel valore critico del PDN oltre il quale l’equilibrio risulta instabile.

• Nel caso riportato in Figura 4.8 si nota chiaramente che per alcuni valori di r1 e r2due autovalori presentano parte reale positiva e di conseguenza l’equilibrio risultainstabile. Si noti inoltre che, a differenza dell’esempio di sezione 4.2.4 ove, sottoopportune condizioni sui parametri, da un certo valore critico in poi di r2 l’equilibrioe sempre instabile, in questo caso sembra esistere una regione limitata di coppie r1ed r2 che rendono l’equilibrio stesso instabile. Tuttavia sembra che per valori dir1 ed r2 sufficientemente piccoli o sufficientemente grandi l’equilibrio sia localmenteasintoticamente stabile.

Invitiamo il lettore a notare come le due situazioni riportate siano esattamente analoghea quelle che proponiamo nell’esempio della sezione 4.2.4 anche se in questo caso non estato possibile determinare le relazioni che consentono di discernere tra i due possibilicasi. Le Figure 4.10 e 4.11 riportano una particolare soluzione del sistema nel caso in cui

Analisi dell’equilibrio (4.2.1d)(coesistenza di entrambe le specie) 63

Figura 4.7: Grafici relativi alla parte reale degli autovalori della matrice Jacobiana (4.2.2) valu-tata nell’equilibrio (4.2.1d) al variare dei due PDN. In particolare notiamo come con la scelta deiparametri effettuata in questo caso gli autovalori siano sempre a parte reale negativa ovvero l’equi-librio (4.2.1d) sia sempre localmenente asintoticamente stabile. Sugli assi orizzontali sono riportatii valori dei PDN r1 e r2 mentre sull’asse verticale la parte reale degli autovalori. I parametri sceltisono: a = 1, b = 0.5, g = 0.2, α = 0.5, β = 1, γ = 1, k1 = 1, k2 = 1


Figura 4.8: Grafici relativi alla parte reale degli autovalori della matrice Jacobiana (4.2.2) valu-tata nell’equilibrio (4.2.1d) al variare dei due PDN. In particolare notiamo come con la scelta deiparametri effettuata in questo caso esistano delle coppie di valori r1 e r2 tali da rendere la partereale di due autovalori positiva. Sugli assi orizzontali sono riportati i valori dei PDN r1 e r2 mentresull’asse verticale la parte reale degli autovalori. I parametri scelti sono: a = 1, b = 0.2, g = 0.8,α = 1, β = 0.5, γ = 1, k1 = 1, k2 = 1. Nella Figura 4.9 e riportata in dettaglio la regione cherende instabile l’equilibrio (4.2.1d).

Un esempio modellisticamente significativo 65

Figura 4.9: Grafico della regione che rende instabile l’equilibrio (4.2.1d) nel piano r1r2. All’internodella regione ci sono le coppie (r1, r2) che rendono l’equilibrio instabile mentre all’esterno quelleche lo rendono localmente asintoticamente stabile. In ascissa sono riportati i valori di r1 mentrein ordinata di r2. I parametri scelti sono: a = 1, b = 0.2, g = 0.8, α = 1, β = 0.5, γ = 1, k1 = 1,k2 = 1 come nei grafici della parte reale degli autovalori riportati in Figura 4.8.

l’equilibrio (4.2.1d) risulti attrattivo. In particolare si sono scelti gli stessi parametri deigrafici relativi agli autovalori associati a tale equilibrio.

4.2.4 Un esempio modellisticamente significativo

Proponiamo qui di seguito un esempio modellisticamente significativo che consente ditrattare l’equilibrio (4.2.1d) (coesistenza di prede e predatori) dal punto di vista analiticomediante le tecniche di linearizzazione. Si consideri il sistema (4.1.1) con i parametri k1, γe β nulli. In tale situazione il sistema prende la forma: x′(t) = x(t)

(a− bx(t)− gy(t)

)y′(t) = y(t)

(− α+ k2

∫ tt0

1r2

exp− t−s

r2

x(s)ds

) .

Si suppone inoltre che tutti gli altri parametri siano strettamente positivi (eccetto il con-trollo logistico b che puo anche pensarsi nullo). Tale sistema risulta molto interessante dalpunto di vista modellistico per i seguenti due motivi:

• Il termine di memoria e contemplato solo nell’interazione tra i predatori e le prede.Tale ipotesi trova riscontro nel fatto che si suppone che i predatori si riproducano


Figura 4.10: Soluzione numerica del sistema (4.1.1) nel caso piu generale atta a mostrare l’at-trattivita dell’equilibrio (4.2.1d). Come funzioni iniziali si sono scelte le costanti x0(t) = 0.5 ey0(t) = 0.5. In particolare in ordinata sono mostrati gli andamenti del numero di prede (rosso) edi predatori (blu) mentre in ascissa il tempo. I parametri relativi a questa simulazione sono: a = 1,b = 0.5, g = 0.2, α = 0.5, k1 = k2 = 1, β = 1, γ = 1. Come si nota tali parametri sono gli stessi deigrafici riportanti la parte reale degli autovalori di Figura 4.7. I parametri r1, r2 sono stati scelti inmodo arbitrario non esistendo, in questo caso, alcuna coppia che renda l’equilibrio instabile. Talescelta e stata r1 = 5 e r2 = 2.

Figura 4.11: Soluzione numerica del sistema (4.1.1) nel caso piu generale atta a mostrare l’at-trattivita dell’equilibrio (4.2.1d). Come funzioni iniziali si sono scelte le costanti x0(t) = 0.5 ey0(t) = 0.5. In particolare in ordinata sono mostrati gli andamenti del numero di prede (rosso) edi predatori (blu) mentre in ascissa il tempo. I parametri relativi a questa simulazione sono: a = 1,b = 0.2, g = 0.8, α = 1, k1 = k2 = 1, β = 0.5, γ = 1. Come si nota tali parametri sono gli stessi deigrafici riportanti la parte reale degli autovalori di Figura 4.8. Si e provveduto inoltre a scegliereuna coppia di valori r1, r2 che non rendesse l’equilibrio instabile in riferimento alla Figura 4.9. Talescelta e stata r1 = 6 e r2 = 0.2.

Un esempio modellisticamente significativo 67

in funzione del numero di prede catturate nel passato mentre che le prede diminui-scano istantaneamente quando vengono catturate dai predatori. Come si nota unamodellizzazione di questo tipo appare ragionevole.

• Il termine di controllo logistico riferito ai predatori e per certi versi superfluo. Sele prede non possono superare una certa soglia (e cio viene tenuto in considerazionemediante il termine di controllo logistico sulle prede stesse), i predatori non possonoacquisire il nutrimento necessario per accrescersi oltre una certa quantita. Pertantoil controllo sui predatori e, in un certo senso, gia contemplato nel controllo delleprede

Ponendo a zero i tre parametri citati la matrice Jacobiana (4.2.2) associata all’equilibrio(4.2.1d) si semplifica notevolmente. In particolare l’equilibrio diventa:

xeq4 = ueq4 =α

k2,

yeq4 = veq4 =ak2 − bα

gk2,

mentre la matrice Jacobiana ad esso associata risulta:

J(α

k2,ak2 − bα

gk2,ak2 − bα

gk2,α

k2) =

− bα

k2−αg

k20 0

0 0 0 ak2−bαg

0 1r1

− 1r1

01r2

0 0 − 1r2

.Il polinomio caratteristico di tale matrice si scrive:

1k2r1r2

(1 + r1λ)(r2k2λ

3 + (k2 + bαr2)λ2 + bαλ+ ak2α− bα2).

Iniziamo col notare che uno degli autovalori della matrice Jacobiana e pari a − 1r1

. Taleautovalore compare, anche se il termine contenente r1 non viene piu considerato, poichela matrice Jacobiana (4.2.2) e associata al sistema nella sua forma piu generale (4.1.3)e non a quella semplificata di questo esempio. Per quest’ultimo sistema si puo infattiintrodurre una sola variabile ausiliaria e di conseguenza andare ad analizzare una matriceJacobiana 3 × 3 in luogo della 4 × 4. Gli autovalori di tali matrici coincidono a meno diλ1 = − 1

r1. Del resto considerare o meno tale autovalore non cambia nulla dal punto di vista

dell’equilibrio in quanto e a parte reale negativa ∀r1 > 0 arbitrariamente fissato. Dal puntodi vista della stabilita dell’equilibrio e invece necessario riuscire a determinare se tutti irestanti tre autovalori hanno parte reale negativa o positiva. Non potendo determinareesplicitamente i tre autovalori λ2, λ3, λ4 sfruttiamo il criterio di Routh-Hurwitz (si vedaappendice A.1.2) per determinare se e quanti autovalori hanno parte reale positiva e inquali casi. Costruiamo la tabella di Routh-Hurwitz associata al polinomio:

r2k2λ3 + (k2 + bαr2)λ2 + bαλ+ ak2α− bα2.

Tale tabella e: r2k2 bα 0

k2 + bαr2 ak2α− bα2 0α(bk2+b2αr2−r2k2

2a+bαr2k2)k2+bαr2

0 0ak2α− bα2 0 0

.

Sulla base del criterio di Routh-Hurwitz il numero di autovalori positivi e pari al numero dicambi di segno nella prima colonna della tabella. I primi due elementi della prima colonnasono sempre strettamente positivi. Il quarto elemento e pure strettamente positivo intutti i casi di interesse modellistico3. Di conseguenza tutti gli autovalori hanno parte realenegativa se:

α(bk2 + b2αr2 − r2k22a+ bαr2k2)

k2 + bαr2> 0,

mentre ve ne sono due con parte reale positiva in caso contrario. In particolare si possonopresentare due diverse situazioni: se b2α+bαk2−k2

2a > 0 al fine di avere tutti gli autovaloricon parte reale negativa deve essere:

r2 > − bk2

b2α+ bαk2 − k22a,

ma cio e sempre verificato in quanto r2 > 0. Se invece si ha: b2α + bαk2 − k22a < 0 la

condizione al fine di avere tutti gli autovalori con parte reale negativa diventa:

0 < r2 < − bk2

b2α+ bαk2 − k22a

e di conseguenza esiste un valore critico del PDN oltre il quale l’equilibrio diventa instabile.Riassumiamo quanto detto nel seguente risultato:

Teorema 4.2.2. Nel caso dell’esempio (k1 = γ = β = 0) allora:

• Se b2α + bαk2 − k22a > 0 l’equilibrio (4.2.1d) (coesistenza di prede e predatori) e

sempre localmente asintoticamente stabile.

• b2α + bαk2 − k22a < 0 l’equilibrio (4.2.1d) (coesistenza di prede e predatori) risulta

localmente asintoticamente stabile se e solo se 0 < r2 < − bk2

b2α+bαk2−k22a

, viceversa einstabile.

Prima di concludere l’esempio analizziamo piu in dettaglio la condizione sotto la qualel’equilibrio risulta sempre localmente asintoticamente stabile. In particolare riscriviamotale condizione in modo tale da poterle conferire un significato modellistico. Ricordiamoche l’equilibrio (4.2.1d), e modellisticamente interessante in questo caso particolare se:

α <a

bk2.

Inoltre e sempre attrattivo se:

b2α+ bαk2 − k22a > 0.

Dividendo l’ultima disequazione per gk22 > 0 e riarrangiando i membri in modo opportuno

si ha:b2

gk2xeq4 > yeq4 .

3Infatti se ak2 − bα ≤ 0 l’equilibrio ha una coordinata negativa o coincide con l’equilibrio (4.2.1b).

Limitatezza delle soluzioni 69

Figura 4.12: Soluzione numerica del sistema (4.1.1) nel caso semplificato proposto nell’esempio(γ = k1 = β = 0) atta a mostrare l’attrattivita dell’equilibrio (4.2.1d). Come funzioni inizialisi sono scelte le costanti x0(t) = 0.5 e y0(t) = 0.5. In particolare in ordinata sono mostrati gliandamenti del numero di prede (rosso) e di predatori (blu) mentre in ascissa il tempo. I parametrirelativi a questa simulazione sono: a = 1, b = 1, g = 0.8, α = 0.9, k2 = 1, r2 = 2. La scelta di r2non influisce sull’attrattivita dell’equilibrio in quanto, in questo caso, non esiste un valore criticodel PDN oltre il quale l’equilibrio risulta instabile.

Cio si puo interpretare dicendo che l’equilibrio di coesistenza tra le specie e model-listicamente interessante se la capacita dell’ambiente per le prede e sufficien-temente grande. Inoltre e sempre attrattivo se i valori di equilibrio prevedonoun numero di prede sufficientemente grande in relazione a quello dei predatori,viceversa e attrattivo solo se la memoria dell’interazione tra specie e corta. Adogni modo si noti che qualora la capacita dell’ambiente sia sufficientemente grande l’equi-librio e sempre localmente asintoticamente stabile se r → 0. E verosimile che nel caso piugenerale, trattabile solo numericamente, si abbiano delle relazioni simili che discrimininoi due casi. Concludiamo l’esempio proponendo le due Figure 4.12 e 4.13 che mostranogli andamenti delle componenti della soluzione nel caso di equilibrio attrattivo sia quandoesiste il valore critico del PDN sia quando tale valore non esiste.

4.3 Limitatezza delle soluzioni

Una proprieta auspicabile per il sistema (4.1.1) e che le soluzioni rimangano limitate neltempo. Tale caratteristica, per quanto di per se non consenta di estendere i risultatidi attrattivita locale a tutto il dominio considerato, assicura perlomeno che le soluzioniesistano per ogni tempo t. Non esistono risultati generali in materia se non legati alla teoriadelle funzioni di Lyapunov. Tuttavia, trovare una funzione di Lyapunov per il sistemaautonomo (4.1.3) che garantisca la limitatezza della sua soluzione e di conseguenza quelladel sistema (4.1.1) non e affatto banale. Un approccio piu diretto, basato sull’analisi dellederivate del sistema integro-differenziale, ci ha tuttavia consentito di formulare i seguentirisultati:

Figura 4.13: Soluzione numerica del sistema (4.1.1) nel caso semplificato proposto nell’esempio(γ = k1 = β = 0) atta a mostrare l’attrattivita dell’equilibrio (4.2.1d). Come funzioni inizialisi sono scelte le costanti x0(t) = 0.5 e y0(t) = 0.5. In particolare in ordinata sono mostrati gliandamenti del numero di prede (rosso) e di predatori (blu) mentre in ascissa il tempo. I parametrirelativi a questa simulazione sono: a = 1, b = 0.3, g = 0.8, α = 1, k2 = 1, r2 = 0.3. La sceltadi r2 influisce sull’attrattivita dell’equilibrio in quanto, in questo caso, esiste un valore critico delPDN oltre il quale l’equilibrio e instabile. In particolare l’equilibrio (4.2.1d) risulta attrattivo solose r2 < 0.4918.

Teorema 4.3.1. Si consideri il sistema di Lotka-Volterra con memoria (4.1.1). Se b > 0la componente x(t) della soluzione e limitata. In particolare se x(0) ≤ a

b allora x(t) ≤ab ∀t > 0 mentre se x(0) = x > a

b si ha che x(t) ≤ x(0) ∀t > 0. Inoltre in questo ultimocaso fissato un arbitrario ε > 0 esiste t > 0 tale che x(t) < a

b + ε ∀t > t.

Dimostrazione. Si consideri x′(t):

x′ = x(a− bx− gy − k1

∫ t

t0

1r1

exp− t− s

r1

y(s)ds

).

Notiamo che per x ≥ ab si ha x′ ≤ 0, potendo l’uguaglianza valere solo se x = a

b :

x′ = x(a− bx− gy − k1

∫ t

t0

1r1

exp− t− s

r1

y(s)ds

)≤ x

(a− bx) ≤ 0.

La dimostrazione si puo a questo punto ritenere conclusa. Infatti se x(0) ≤ ab la soluzione

non puo mai superare tale valore (infatti per superarlo dovrebbe avere derivata positivain a

b ). Viceversa se x(0) > ab la derivata della soluzione si mantiene negativa almeno

finche non si ha x = ab . Di conseguenza x(t) e decrescente e fissato un arbitrario valore

ε > 0 esiste sicuramente un opportuno istante temporale t dipendente da ε tale che x(t) <ab + ε ∀t > t.

Teorema 4.3.2. Si consideri il sistema di Lotka-Volterra con memoria (4.1.1). Sianob > 0 e β > 0. Allora la soluzione e limitata su tutto il semiasse temporale positivo.

Limitatezza delle soluzioni 71

Dimostrazione. Grazie al Teorema 4.3.1 la componente x(t) della soluzione e limitata neltempo e in particolare se x(0) ≤ a

b allora x(t) ≤ ab , viceversa esistono comunque un ε > 0

e un istante t tale che x(t) < ab + ε ∀t > t. Tenuto conto di cio la derivata della y si puo

maggiorare da t in poi nel seguente modo:

y′ = y(− α− βy + γx+ k2

∫ t

t0

1r2

exp− t− s

r2

x(s)ds

)≤ y(− α− βy + γ(

a

b+ ε) + k2

∫ t

t0

1r2

exp− t− s

r2

(a

b+ ε)ds

)≤ y(−αb+ γa+ k2a− βy + (γ + k2)ε

b

),

e l’ultimo membro e negativo ogniqualvolta si ha:

y >−αb+ γa+ k2a+ (γ + k2)ε

bβ.

Data l’arbitrarieta di ε si puo concludere che se y(0) e minore di −αb+γa+k2abβ la componente

y(t) non supera mai tale valore, mentre se y(0) e maggiore −αb+γa+k2abβ valore critico la

componente della soluzione tende a riportarsi sul o sotto tale valore e cio conclude.

Teorema 4.3.3. Si consideri il sistema di Lotka-Volterra con memoria (4.1.1). Se b > 0e l’equilibrio (4.2.1b) e attrattivo allora tale equilibrio e anche globalmente asintoticamentestabile.

Dimostrazione. Iniziamo col ricordare che se l’equilibrio (4.2.1b) e attrattivo allora k2a+γa − αb < 0. Grazie al Teorema 4.3.1 la componente x(t) della soluzione e limitata neltempo e in particolare se x(0) ≥ a

b allora x(t) ≤ ab , viceversa esistono comunque un ε > 0

e un istante t tale che x(t) < ab + ε ∀t > t. Tenuto conto di cio la derivata della y si puo

maggiorare nel seguente modo:

y′ = y(− α− βy + γx+ k2

∫ t

t0

1r2

exp− t− s

r2

x(s)ds

)≤ y(− α+ γ(

a

b+ ε) + k2

∫ t

t0

1r2

exp− t− s

r2

(a

b+ ε)ds

)≤ y(k2a+ γa− αb

b+ (γ + k2)ε

).

Poiche k2a+γa−αbb e strettamente negativa, e sempre possibile scegliere ε sufficientemente

piccolo tale che da t in poi si abbia:

k2a+ γa− αb

b+ (γ + k2)ε < 0.

Si noti inoltre che la componente y(t) della soluzione non puo esplodere prima di t grazie

alla limitatezza della x. Infatti y′ si mantiene limitata nel tempo:

y′ = y(− α− βy + γx+ k2

∫ t

t0

1r2

exp− t− s

r2

x(s)ds

)≤ y(− α+ γx+ k2

∫ t

t0

1r2

exp− t− s

r2

xds)

≤ y((γ + k2)x

)< +∞,

ove x = maxx(t). Poiche da t in poi si ha y′ < 0 si ha che la componente y della soluzionetende a zero. In particolare fissato ε1 > 0 arbitrario, esiste t1 tale che y < ε1 ∀t > t1. Maquesto implica che sia sempre possibile scegliere ε1 tale che, fissato un ε2 > 0 arbitrario,da t1 in poi valga:

gy + k1

∫ t

t0

1r1

exp− t− s

r1

y(s) < ε2.

Ma allora per t > t1 e x < a−ε2b si ha:

x′ = x(a− bx− gy − k1

∫ t

t0

1r1

exp− t− s

r1

y(s)ds

)> x

(a− bx− ε2) ≥ 0,

e per l’arbitrarieta di ε e di ε2 si conclude che anche la x deve tendere al valore di equilibrioab .

Gli enunciati dei Teoremi appena proposti si possono apprezzare in relazione alle simula-zioni numeriche gia presentate. In particolare per il Teorema 4.3.1 si veda una qualunqueFigura corrispondente a una simulazione con b > 0, per il Teorema 4.3.3 si vedano leFigure 4.3 e 4.5 ed infine, per il Teorema 4.3.2 si vedano le Figure 4.10 e 4.11. I risultatinumerici proposti e quelli ottenuti durante lo svolgimento del nostro lavoro farebbero sup-porre che valga un risultato piu forte per quanto riguarda la limitatezza delle soluzioni.In particolare sembra che anche nel caso piu generale non sia necessario che β > 0 pergarantire che entrambe le componenti siano limitate ma che basti considerare positivo ilparametro di controllo logistico sulle prede. Tale ipotesi e anche supportata dal ragiona-mento di carattere modellistico gia proposto nell’esempio. Infatti i predatori aumentanoin funzione delle prede assimilate mentre le prede diminuiscono in funzione del numero dipredatori. In particolare se le prede sono limitate i predatori non potranno assimilare cibosufficiente per crescere in modo illimitato. Questo ragionamento euristico, affiancato alleprove numeriche citate, ci ha spinto a formulare la congettura seguente:

Congettura 4.3.1. Qualora b > 0 le componenti della soluzione si mantengono limitatenel tempo.

In ogni caso nell’ipotesi b > 0 il numero di prede non puo mai superare una certa soglia(grazie al Teorema 4.3.1) e il numero di predatori non puo divergere in modo monotono.Questa ultima affermazione e giustificata dal seguente Teorema:

Cicli limite 73

Teorema 4.3.4. Qualora b > 0 la variabile y non puo divergere in modo definitivamentemonotono.

Dimostrazione. Iniziamo col ricordare che in base al Teorema 4.3.1 comunque si scelgax(0) e un ε > 0 esiste un t tale che: x ≤ a

b + ε per ogni t > t. Inoltre se β > 0 il Teorema4.3.2 garantisce che entrambe le soluzioni siano limitate e di conseguenza l’unico caso incui la componente y della soluzione puo divergere e quello in cui β = 0. Supponiamo perassurdo che la soluzione y diverga in modo definitivamente monotono da un certo istantet1 in poi. Cio implica che esiste un certo istante t2 tale che y > a−bx

g+k1per ogni t > t2.

Tuttavia se y e maggiore di a−bxg+k1

si ha x′ < 0. Infatti se y = maxt≤t2 y(t) si ha:

x′ = x(a− bx− gy − k1

∫ t

t0

1r1

exp− t− s

r1

y(s)ds

)≤ x

(a− bx− gy − k1

∫ t

t0

1r1

exp− t− s

r1

yds)

≤ x(a− bx− gy − k1y) < 0

Ma allora, essendo y > a−bxg+k1

per ogni t > t2 si ha anche x′ < 0 per ogni t > t2. Tenendoconto che il nucleo scelto e decrescente e pertanto tiene traccia prevalentemente dei valoripiu recenti della variabile esiste un certo istante t3 dopo il quale si ha:

γx+ k2

∫ t

t0

1r1

exp− t− s

r2

x(s)ds < α

Ma cio implica y′ < 0, assurdo.

4.4 Cicli limite

Nelle sezioni 4.2 e 4.3 abbiamo illustrato sotto quali condizioni gli equilibri del sistema(4.1.3) (e i corrispettivi equilibri del sistema (4.1.1)) sono attrattivi o repulsivi e sottoquali condizioni le soluzioni si mantengono limitate. Tralasciando il caso con β = 0per cui la limitatezza della soluzione e stata solo congetturata ma non dimostrata, cisoffermiamo ora sul caso in cui β > 0. In tale caso, supposto ovviamente b > 0, ilTeorema 4.3.2 garantisce che la soluzione sia limitata nel tempo. Tuttavia, nell’analisidegli equilibri fatta nella sezione 4.2 abbiamo notato come possa succedere che tutti gliequilibri siano repulsivi. In particolare tale situazione si verifica quando l’equilibrio (4.2.1d)ha entrambe le coordinate positive e si scelgono i parametri in modo tale che esista unacoppia di valori r1 ed r2 che rendono instabile l’equilibrio medesimo. Ci si chiede qualesia il comportamento della soluzione in questo caso particolare. Le ipotesi che si possonoformulare sono essenzialmente due:

• la soluzione ha un comportamento caotico,

• la soluzione si stabilizza su un ciclo limite.


Figura 4.14: Soluzione numerica atta a mostrare lo svilupparsi di un ciclo limite nel caso in cuile componenti della soluzione siano limitate e non esista alcun equilibrio attrattivo. In particolarel’equilibrio (4.2.1d) risulta instabile a causa della particolare scelta di r1 ed r2. Il grafico riportain ordinata gli andamenti in funzione del tempo del numero di prede (rosso) e di predatori (blu)mentre in ascissa il tempo. I parametri scelti sono: a = 1, b = 0.2, g = 0.8, α = 1, β = 0.5,γ = 1, k1 = 1, k2 = 1 come in Figura 4.8. Si noti in Figura 4.9 come la particolare scelta effettuatar1 = r2 = 1 renda l’equilibrio instabile. Come funzioni iniziali si e scelta la coppia di costantix0(t) = y0(t) = 0.5.

Le prove numeriche fanno propendere per la seconda ipotesi. In particolare formalizziamotale supposizione nella seguente Congettura:

Congettura 4.4.1. Si consideri il sistema di Lotka-Volterra con memoria (4.1.1). Se lesoluzioni sono limitate e tutti gli equilibri di interesse modellistico sono repulsivi alloraogni traiettoria del sistema con dati iniziali positivi tende a stabilizzarsi su un ciclo limite.

In questo senso proponiamo alcune simulazioni numeriche relative al caso in cui esistauna coppia di valori di PDN che rendono instabile l’equilibrio (4.2.1d). Le Figure 4.14,4.15, 4.16 e 4.17 illustrano come si sviluppi un ciclo limite nel caso l’equilibrio sia instabilecon i parametri relativi al grafico degli autovalori della matrice Jacobiana dell’equilibrio(4.2.1d) mostrati in Figura 4.8. Tale ciclo limite sembra svilupparsi sia nel caso in cuile coordinate iniziali siano all’interno del ciclo stesso sia nel caso in cui siano all’esterno.Come anticipato nelle motivazioni che ci hanno spinto a formulare la Congettura 4.3.1,dalle simulazioni numeriche sembra che le soluzioni siano limitate anche qualora β = 0.In particolare anche in tale caso i risultati ottenuti numericamente sembrano confermarela possibilita che si sviluppi un ciclo limite. In particolare nelle Figure 4.18 e 4.19 vieneriportata una soluzione numerica del sistema semplificato proposto nell’esempio di sezione4.2.4 nel caso in cui il PDN renda l’equilibrio instabile. Si noti che, in linea di principio,esiste un ciclo limite ma non esiste alcuna orbita periodica. Infatti se vi fosse un’ orbitaperiodica per il sistema 4.1.1 dovrebbe esistere una orbita periodica anche per il sistema4.1.3. Ma cio almeno in generale non si verifica. Infatti abbiamo mostrato nella sezione 2.5come la convoluzione tra il nucleo (2.2.1) e una funzione periodica sembra, a partire dallesimulazioni numeriche, assestarsi su un ciclo periodico ma non e, in generale, periodica.

Cicli limite 75

Figura 4.15: Soluzione numerica atta a mostrare lo svilupparsi di un ciclo limite nel caso in cuile componenti della soluzione siano limitate e non esista alcun equilibrio attrattivo. In particolarel’equilibrio (4.2.1d) risulta instabile a causa della particolare scelta di r1 ed r2. Il grafico riportala traiettoria della soluzione nel piano xy (in ordinata i predatori mentre in ascissa le prede). Iparametri scelti sono: a = 1, b = 0.2, g = 0.8, α = 1, β = 0.5, γ = 1, k1 = 1, k2 = 1 come inFigura 4.8. Si noti in Figura 4.9 come la particolare scelta effettuata r1 = r2 = 1 renda l’equilibrioinstabile. Come funzioni iniziali si e scelta la coppia di costanti x0(t) = y0(t) = 0.5.

Figura 4.16: Soluzione numerica atta a mostrare lo svilupparsi di un ciclo limite nel caso in cuile componenti della soluzione siano limitate e non esista alcun equilibrio attrattivo. In particolarel’equilibrio (4.2.1d) risulta instabile a causa della particolare scelta di r1 ed r2. Il grafico riportain ordinata gli andamenti in funzione del tempo del numero di prede (rosso) e di predatori (blu)mentre in ascissa il tempo. I parametri scelti sono: a = 1, b = 0.2, g = 0.8, α = 1, β = 0.5,γ = 1, k1 = 1, k2 = 1 come in Figura 4.8. Si noti in Figura 4.9 come la particolare scelta effettuatar1 = r2 = 1 renda l’equilibrio instabile. Come funzioni iniziali si e scelta la coppia di costantix0(t) = y0(t) = 5.


Figura 4.17: Soluzione numerica atta a mostrare lo svilupparsi di un ciclo limite nel caso in cuile componenti della soluzione siano limitate e non esista alcun equilibrio attrattivo. In particolarel’equilibrio (4.2.1d) risulta instabile a causa della particolare scelta di r1 ed r2. Il grafico riportala traiettoria della soluzione nel piano xy (in ordinata i predatori mentre in ascissa le prede). Iparametri scelti sono: a = 1, b = 0.2, g = 0.8, α = 1, β = 0.5, γ = 1, k1 = 1, k2 = 1 come inFigura 4.8. Si noti in Figura 4.9 come la particolare scelta effettuata r1 = r2 = 1 renda l’equilibrioinstabile. Come funzioni iniziali si e scelta la coppia di costanti x0(t) = y0(t) = 5.

Figura 4.18: Soluzione numerica atta a mostrare lo svilupparsi di un ciclo limite nel caso sem-plificato proposto nell’esempio 4.2.4. In particolare l’equilibrio (4.2.1d) risulta instabile a causadella particolare scelta r2 = 0.55 > 0.4918 = rcritico. Il grafico riporta in ordinata gli andamenti infunzione del tempo del numero di prede (rosso) e di predatori (blu) mentre in ascissa il tempo. Iparametri scelti sono: a = 1, b = 0.3, g = 0.8, α = 1, β = 0, γ = 0, k1 = 0, k2 = 1. Come funzioniiniziali si sono scelte le costanti x0(t) = 0.2 e y0(t) = 0.4.

Cicli limite 77

Figura 4.19: Soluzione numerica atta a mostrare lo svilupparsi di un ciclo limite nel caso sem-plificato proposto nell’esempio 4.2.4. In particolare l’equilibrio (4.2.1d) risulta instabile a causadella particolare scelta r2 = 0.55 > 0.4918 = rcritico. Il grafico riporta la traiettoria della soluzionenel piano xy (in ordinata i predatori mentre in ascissa le prede). I parametri scelti sono: a = 1,b = 0.3, g = 0.8, α = 1, β = 0, γ = 0, k1 = 0, k2 = 1. Come funzioni iniziali si sono scelte lecostanti x0(t) = 0.2 e y0(t) = 0.4.

Nel caso in cui si presenta un ciclo limite, integrando numericamente le soluzioni del sistemasu un “periodo”4 una volta che l’orbita si sia sufficientemente stabilizzata, sembrerebbe chela media integrale su tale “periodo” di ciascuna delle due variabili eguagli il corrispettivovalore dell’equilibrio (4.2.1d). Enunciamo ora un risultato, da noi ottenuto, che dimostratale supposizione.

Teorema 4.4.1. Si consideri il sistema (4.1.1) e si supponga che la soluzione tenda a unciclo limite. Sia t > 0 un opportuno tempo positivo oltre il quale si puo pensare che lasoluzione sia diventata periodica di periodo T . Allora per qualunque tempo t > t valgonole seguenti due relazioni:

x(t) =1T

∫ t+T

tx(s)ds =

aβ + αk1 + αg

gk2 + γk1 + gγ + bβ + k1k2+ ou,v,y,x(t), (4.4.1a)

y(t) =1T

∫ t+T

ty(s)ds =

aγ + ak2 − bα

gk2 + γk1 + gγ + bβ + k1k2+ ou,v,y,x(t). (4.4.1b)

ove ou,v,y,x(t) e una quantita che tende a zero quando t→ +∞ dipendente dalla particolaresoluzione analizzata.

Dimostrazione. Si consideri il sistema (4.1.3) associato al sistema (4.1.1) e, supposto chetutte e quattro le variabili si siano sufficientemente stabilizzate su un ciclo periodico siintegrino tutte le equazioni sul periodo:∫ t+T

t

x′(s)x(s)

ds = aT − b

∫ t+T

tx(s)ds− g

∫ t+T

ty(s)ds− k1

∫ t+T

tv(s)ds,

4In realta poiche l’orbita non e periodica ma tende ad un’orbita periodica l’intervallo temporale diintegrazione tra, per esempio, due picchi consecutivi di una delle due variabili non puo ritenersi strettamenteun periodo. Con abuso di linguaggio utilizzeremo ugualmente tale termine nel seguito di questa sezione.


∫ t+T

t

y′(s)y(s)

ds = −αT − β

∫ t+T

ty(s)ds+ γ

∫ t+T

tx(s)ds+ k2

∫ t+T

tu(s)ds,∫ t+T

tv′(s)ds =

1r1

(∫ t+T

ty(s)ds−

∫ t+T

tv(s)ds

),∫ t+T

tu′(s)ds =

1r2

(∫ t+T

tx(s)ds−

∫ t+T

tu(s)ds

).

Dalle ultime due relazioni, considerando che u(t+T ) = u(t)+ou(t) e v(t+T ) = v(t)+ov(t)si ha: ∫ t+T

ty(s)ds =

∫ t+T

tv(s)ds+ ov(t),∫ t+T

tx(s)ds =

∫ t+T

tu(s)ds+ ou(t).

Tenendo inoltre conto che∫ t+Tt

y′(s)y(s) ds = oy(t) e

∫ t+Tt

x′(s)x(s) ds = ox(t) si puo arrivare a

scrivere il seguente sistema di due equazioni in due incognite per x(t) e y(t):a− bx− gy − k1y = ov,x(t)−α− βy + γx+ k2x = ou,y(t)

,

che risolto fornisce il risultato voluto.

Concludiamo questa sezione sui cicli limite notando l’estrema rilevanza che un compor-tamento di tale tipo ha dal punto di vista modellistico. Il grande pregio del sistema diLotka-Volterra classico e quello di avere orbite periodiche che spesso possono modelliz-zare in modo appropriato i reali andamenti del numero di individui di una popolazionenel caso di interazione preda-predatore. Tuttavia tale modello presenta delle limitatezzelegate alla sua semplicita. L’ipotesi di interazione istantanea preda-predatore e di assenzadi controllo logistico (ovvero di risorse illimitate) sono evidentemente semplificazioni nonaccettabili. La possibilita di riscontrare andamenti di tipo periodico in un sistema conte-nente un’interazione non piu solo istantanea ma anche contenente termini di “memoria”nonche i controlli logistici appare pertanto estremamente significativa. Questa possibilitasembra poter sussistere nel caso da noi analizzato e conferisce pertanto a tale caso unagrande importanza.

4.5 Soluzioni al variare dei PDN

In questa ultima sezione ci occupiamo di trattare come le soluzioni variano al variare deidue PDN. In particolare ci concentriamo su due casi notevoli: quello in cui entrambi iPDN tendano a zero e quello in cui entrambi i PDN tendano a infinito.

4.5.1 Convergenza del sistema a quello senza memoria quando il PDNtende a zero

Questo caso e quello di interesse modellistico piu importante e coincide con la situazionein cui vi e perdita di memoria: la storia delle variabili e influenzata solo dal valore corrente

Soluzioni al variare dei PDN 79

delle medesime e il sistema, almeno formalmente degenera in un sistema di Lotka-Volterracon controllo logistico con opportuni parametri. Infatti si ricordi che il nucleo (2.2.1)converge nel senso delle distribuzioni alla delta di Dirac. Sostituendo il nucleo con la deltanel sistema (4.1.1) si ha: x′(t) = x(t)

(a− bx(t)− gy(t)− k1

∫ tt0δ(t− s)y(s)ds

)y′(t) = y(t)

(− α− βy(t) + γx(t) + k2

∫ tt0δ(t− s)x(s)ds

) ,

che formalmente si riduce a: x′(t) = x(t)(a− bx(t)− (g + k1)y(t)

)y′(t) = y(t)

(− α− βy(t) + (γ + k2)x(t)

) . (4.5.1)

Evidentemente la (4.5.1) puo essere vista come un particolare caso del sistema (1.5.1). Adogni modo non vi e alcuna garanzia che le soluzioni del modello di Lotka-Volterra conmemoria tendano a quelle del modello degenere introdotto. Stime sulla convergenza dellesoluzioni a quelle del modello degenere si possono effettuare mediante tecniche che esulanodagli scopi della nostra trattazione. Siamo comunque interessati a una verifica numericadella convergenza uniforme tra le soluzioni della (4.5.1) e quelle del corrispondente sistemadi tipo (1.5.1). Come peraltro gia notato la convergenza uniforme e una convergenza moltoforte che implica quelle in norma integrale. Mostriamo in questo senso le Figure 4.20,4.21, 4.22 e 4.23. Si noti in particolare come le soluzioni gia con valori di r1 ed r2 nonparticolarmente piccoli siano molto simili. Diminuendo ulteriormente i PDN non si riescenemmeno a distinguere i due grafici e pertanto simulazioni relative a valori inferiori deiPDN non vengono riportate.

4.5.2 Convergenza del sistema a quello senza memoria quando il PDNtende a infinito

In questo caso grazie al Microteorema 2.4.2, almeno per tempi limitati, il sistema: x′(t) = x(t)(a− bx(t)− gy(t)− k1

∫ tt0δ(t− s)y(s)ds

)y′(t) = y(t)

(− α− βy(t) + γx(t) + k2

∫ tt0δ(t− s)x(s)ds

)si riconduce formalmente a: x′(t) = x(t)

(a− bx(t)− gy(t)

)y′(t) = y(t)

(− α− βy(t) + γx(t)

) . (4.5.2)

Ad ogni modo non vi e alcuna garanzia che le soluzioni del modello di Lotka-Volterracon memoria tendano a quelle del modello degenere introdotto. Analogamente a quantofatto in precedenza proponiamo una verifica numerica della convergenza uniforme tra lesoluzioni della (4.5.2) e quelle del corrispondente sistema di tipo (1.5.1). Mostriamo inquesto senso le Figure 4.24, 4.25, 4.26 e 4.27. Si noti come una volta fissato un tempofinale T > 0 finito e un ε arbitrario sufficientemente piccolo, sembra esistere sempre unacoppia di valori r1 e r2 tale che le due soluzioni differiscano al piu di ε per ogni t < T .


Figura 4.20: Soluzione numerica atta a mostrare la convergenza della soluzione del sistema (4.1.1)a quella del corrispondente sistema (4.5.1) quando entrambi i PDN tendono a zero. In questo casola soluzione tende all’unico equilibrio attrattivo del sistema (4.5.1). Il grafico riporta in ordinatagli andamenti in funzione del tempo del numero di prede (rosso per il sistema (4.1.1) e giallo peril sistema (4.5.1)) e di predatori (blu per il sistema (4.1.1) e azzurro per il sistema (4.5.1)) mentrein ascissa il tempo. I parametri scelti sono: a = 1, b = 0.2, g = 0.8, α = 1, β = 0.5, γ = 1, k1 = 1,k2 = 1, r1 = r2 = 0.05. Come funzioni iniziali si e scelta la coppia di costanti x0(t) = y0(t) = 0.5.

Figura 4.21: Soluzione numerica atta a mostrare la convergenza della soluzione del sistema (4.1.1)a quella del corrispondente sistema (4.5.1) quando entrambi i PDN tendono a zero. In questo casola soluzione tende all’unico equilibrio attrattivo del sistema (4.5.1). Il grafico riporta le traiettoriedelle soluzioni nel piano xy. In particolare la soluzione del sistema (4.1.1) e riportata in verdementre quella del sistema (4.5.1) in arancio. In ordinata vi e il numero di predatori mentre inascissa il numero di prede. I parametri scelti sono: a = 1, b = 0.2, g = 0.8, α = 1, β = 0.5,γ = 1, k1 = 1, k2 = 1, r1 = r2 = 0.05. Come funzioni iniziali si e scelta la coppia di costantix0(t) = y0(t) = 0.5.


Figura 4.22: Soluzione numerica atta a mostrare la convergenza della soluzione del sistema (4.1.1)a quella del corrispondente sistema (4.5.1) quando entrambi i PDN tendono a zero. In questo casoessendo ambedue i parametri di controllo logistico nulli la soluzione tende a quella del corrispettivosistema di Lotka-Volterra classico. Il grafico riporta in ordinata gli andamenti in funzione deltempo del numero di prede (rosso per il sistema (4.1.1) e giallo per il sistema (4.5.1)) e di predatori(blu per il sistema (4.1.1) e azzurro per il sistema (4.5.1)) mentre in ascissa il tempo. I parametriscelti sono: a = 1, b = 0, g = 0.8, α = 1, β = 0, γ = 1, k1 = 1, k2 = 1, r1 = r2 = 0.05. Comefunzioni iniziali si e scelta la coppia di costanti x0(t) = y0(t) = 0.5.

Figura 4.23: Soluzione numerica atta a mostrare la convergenza della soluzione del sistema (4.1.1)a quella del corrispondente sistema (4.5.1) quando entrambi i PDN tendono a zero. In questo casoessendo ambedue i parametri di controllo logistico nulli la soluzione tende a quella del corrispettivosistema di Lotka-Volterra classico. Il grafico riporta le traiettorie delle soluzioni nel piano xy. Inparticolare la soluzione del sistema (4.1.1) e riportata in verde mentre quella del sistema (4.5.1) inarancio. In ordinata vi e il numero di predatori mentre in ascissa il numero di prede. I parametriscelti sono: a = 1, b = 0, g = 0.8, α = 1, β = 0, γ = 1, k1 = 1, k2 = 1, r1 = r2 = 0.05. Comefunzioni iniziali si e scelta la coppia di costanti x0(t) = y0(t) = 0.5.


Figura 4.24: Soluzione numerica atta a mostrare la convergenza della soluzione del sistema (4.1.1)a quella del corrispondente sistema senza memoria quando entrambi i PDN tendono a infinito. Inquesto caso la soluzione tende all’unico equilibrio attrattivo del sistema senza memoria. Il graficoriporta in ordinata gli andamenti in funzione del tempo del numero di prede (rosso per il sistema(4.1.1) e giallo per il sistema senza memoria) e di predatori (blu per il sistema (4.1.1) e azzurroper il sistema senza memoria) mentre in ascissa il tempo. I parametri scelti sono: a = 1, b = 0.2,g = 0.8, α = 1, β = 0.5, γ = 1, k1 = 1, k2 = 1, r1 = r2 = 10000. Come funzioni iniziali si e sceltala coppia di costanti x0(t) = y0(t) = 0.5.

Figura 4.25: Soluzione numerica atta a mostrare la convergenza della soluzione del sistema (4.1.1)a quella del corrispondente sistema senza memoria quando entrambi i PDN tendono a infinito. Inquesto caso la soluzione tende all’unico equilibrio attrattivo del sistema senza memoria. Il graficoriporta le traiettorie delle soluzioni nel piano xy. In particolare la soluzione del sistema (4.1.1) eriportata in verde mentre quella del sistema senza memoria in arancio. In ordinata vi e il numerodi predatori mentre in ascissa il numero di prede. I parametri scelti sono: a = 1, b = 0.2, g = 0.8,α = 1, β = 0.5, γ = 1, k1 = 1, k2 = 1, r1 = r2 = 10000. Come funzioni iniziali si e scelta la coppiadi costanti x0(t) = y0(t) = 0.5.


Figura 4.26: Soluzione numerica atta a mostrare la convergenza della soluzione del sistema (4.1.1)a quella del corrispondente sistema senza memoria quando entrambi i PDN tendono a infinito. Inquesto caso essendo ambedue i parametri di controllo logistico nulli la soluzione tende a quelladel corrispettivo sistema di Lotka-Volterra classico. Il grafico riporta in ordinata gli andamenti infunzione del tempo del numero di prede (rosso per il sistema (4.1.1) e giallo per il sistema senzamemoria) e di predatori (blu per il sistema (4.1.1) e azzurro per il sistema senza memoria) mentrein ascissa il tempo. I parametri scelti sono: a = 1, b = 0, g = 0.8, α = 1, β = 0, γ = 1, k1 = 1,k2 = 1, r1 = r2 = 10000. Come funzioni iniziali si e scelta la coppia di costanti x0(t) = y0(t) = 0.5.

Figura 4.27: Soluzione numerica atta a mostrare la convergenza della soluzione del sistema (4.1.1)a quella del corrispondente sistema senza memoria quando entrambi i PDN tendono a infinito. Inquesto caso essendo ambedue i parametri di controllo logistico nulli la soluzione tende a quelladel corrispettivo sistema di Lotka-Volterra classico. Il grafico riporta le traiettorie delle soluzioninel piano xy. In particolare la soluzione del sistema (4.1.1) e riportata in verde mentre quella delsistema senza memoria in arancio. In ordinata vi e il numero di predatori mentre in ascissa ilnumero di prede. I parametri scelti sono: a = 1, b = 0, g = 0.8, α = 1, β = 0, γ = 1, k1 = 1,k2 = 1, r1 = r2 = 10000. Come funzioni iniziali si e scelta la coppia di costanti x0(t) = y0(t) = 0.5.

Conclusioni

Vorremmo concludere il nostro lavoro riassumendo schematicamente i maggiori risul-tati che sono stati proposti e ottenuti nell’elaborato. In particolare abbiamo deciso diraggruppare tali risultati per capitolo:

Capitolo 1 In questo capitolo ci siamo occupati di passare in rassegna i modelli biologicifondamentali, in particolare quelli di Malthus, di controllo logistico, di Verhulst-Volterra, di Lotka-Volterra, di Lotka-Volterra con controllo logistico e di Lotka-Volterra con memoria. E un capitolo essenzialmente descrittivo che riassume per-lopiu risultati gia noti. La parte, a nostro avviso, di maggior interesse e quella deiTeoremi e Congetture relativi al sistema di Lotka-Volterra classico contenuti nellepagine 9-14. Questi risultati da noi proposti possono avere applicazioni in camponumerico in quanto, una volta nota una soluzione numerica, consentono di genera-re le soluzioni di sistemi che hanno certe relazioni con quello risolto semplicementemediante dilatazioni. Inoltre la congettura 1.4.1 fornisce una caratterizzazione delperiodo della soluzione in termini di medie integrali delle variabili.

Capitolo 2 In questo capitolo abbiamo introdotto i nuclei di memoria con particolare rife-rimento a quello a decadimento esponenziale (2.2.1). I risultati di maggior interessedi questo capitolo sono relativi alle buone proprieta del nucleo stesso. In particolaresi segnala la convergenza alla delta di Dirac nel senso delle distribuzioni quando ilparametro PDN tende a zero (importante nel considerare sistemi con sempre meno“memoria”) e la possibilita di passare da sistemi integro-differenziali a sistemi solodifferenziali in piu equazioni ma autonomi.

Capitolo 3 In questo capitolo abbiamo trattato in dettaglio l’equazione integro-differenzialedi Verhulst-Volterra nel caso in cui si scelga un nucleo convolutivo del tipo (2.2.1).Risultati importanti riguardano l’assenza di soluzioni periodiche e la attrattivita glo-bale dell’unico equilibrio attrattivo. Altri risultati riguardano la limitatezza dellesoluzioni con speciale riferimento a un risultato proposto da Burton e l’assenza dipossibilita di estinzione qualora il numero iniziale di individui sia positivo. Infinesono stati proposti risultati numerici di convergenza della soluzione dell’equazione diVerhulst-Volterra a quella del modello di controllo logistico qualora il PDN tenda azero o a infinito.

Capitolo 4 In questo capitolo abbiamo trattato in dettaglio il sistema integro-differenzialedi Lotka-Volterra con memoria nel caso in cui si scelga un nucleo convolutivo del tipo

85

86 Conclusioni

(2.2.1). E a nostro avviso il capitolo piu interessante di tutto l’elaborato. Risulta-ti importanti riguardano la limitatezza delle soluzioni, la globale attrattivita di unequilibrio (se sussistono certe relazioni tra i parametri) e la verifica numerica dellapossibile presenza di cicli limite. Questa ultima possibilita e quella che conferisceinteresse a questo capitolo sia per la rilevanza che una tale situazione comporta dalpunto di vista modellistico sia per l’assenza di risultati teorici in questo ambito. Seb-bene non siamo riusciti a dimostrare l’esistenza di un ciclo-limite a livello teoricoriteniamo che il nostro lavoro possa essere interessante e uno spunto per successiviapprofondimenti in materia. Sono stati inoltre fatti dei raffronti tra il modello conmemoria e quelli senza memoria a cui il sistema di Lotka-Volterra con memoria siriconduce formalmente quando i PDN tendono a zero o a infinito. Esempi numericiintegrano tale trattazione verificando che la convergenza non sussiste dal solo pun-to di vista formale ma che anche le soluzioni tendono a quelle dei modelli degenerialmeno per tempi finiti.

Matteo LesinigoAnna Zanzottera

Appendice A

Risultati utilizzati nel testo

In questa appendice vengono riportati i principali risultati che sono stati utilizzati neltesto.

A.1 Algebra Lineare

Teorema A.1.1. Si consideri una matrice quadrata A ∈ R2x2. Allora condizione suffi-ciente affinche tale matrice abbia autovalori a parte reale negativa e che valgano contem-poraneamente le tre disuguaglianze a11 < 0, Tr(A) < 0 e a12 a21 ≤ 0.

Dimostrazione. La dimostrazione e piuttosto semplice. Gli autovalori della matrice Asono:

λ1 =12(a11 + a22 −

√(a11 − a22)2 + 4 a12 a21)

λ2 =12(a11 + a22 +

√(a11 − a22)2 + 4 a12 a21)

Evidentemente se Tr(A) < 0 si ha che la parte reale di λ1 e essa stessa negativa in quantoRe(λ1) ≤ Re(a11 +a22) = Re(Tr(A)). Per quanto riguarda λ2 il discorso e leggermente piucomplesso. Tuttavia se a12 a21 ≤ 0 come ipotizzato si puo effettuare la seguente catena dimaggiorazioni:

λ2 =12(a11 + a22 +

√(a11 − a22)2 + 4 a12 a21)

≤ 12(a11 + a22 +

√(a11 − a22)2)

≤ 12(a11 + a22 + (a11 − a22)) = a11 < 0

Teorema A.1.2 (Criterio di Routh-Hurwitz). Sia:

∆A(λ) = a0λn + a1λ

n−1 + ...+ an

87

88 Appendice

il polinomio caratteristico associato a una matrice A. Si costruisca la seguente tabella(detta tabella di Routh) di dimensioni (n+ 1)× (n+ 1)

r00 r01 r02 . . . r0n

r10 r11 r12 . . . r1n

r20 r21 r22 . . . r2n...

......

. . ....

rn0 rn1 rn2 . . . rnn

dove gli elementi della prima e della seconda riga sono, rispettivamente, i coefficienti conindice pari (a0, a2,...) e dispari (a1, a3,...) del polinomio caratteristico e ai = 0 per i > nmentre tutti gli altri elementi vanno calcolati con la formula:

ri+1,j = − 1ri0

det[ri−1,0 ri−1,j+1

ri,0 ri,j+1

].

Condizione necessaria e sufficiente affinche tutti gli autovalori di A abbiano parte realenegativa e che non vi siano cambi di segno ne elementi nulli nella prima colonna dellatabella di Routh. Inoltre se nella prima colonna non vi sono elementi nulli il numero deglieventuali cambi di segno coincide con il numero degli autovalori a parte reale positiva (Cfr:[DPRFA, pag. 44]).

A.2 Analisi Funzionale

Teorema A.2.1 (Lebesgue). Sia un una successione di funzioni complesse integrabiliconvergente q.o. ad una funzione u definita q.o. in Rn. Se esiste una funzione φ integrabiletale che:

|un(x)| ≤ φ(x) q.o.∀n,

allora u e integrabile e valgono le formule:

limk→∞

∫uk =

∫u

limk→∞

∫|uk − u| = 0

Dimostrazione. Per la dimostrazione si veda [GA3, pagg. 43 e 52]

Definizione A.2.1 (Convergenza in D′). Siano uk una successione di distribuzioni inΩ e u una distribuzione in Ω. Si dice che uk converge ad u nel senso delle distribuzioniuk → u in D′(Ω) quando:

limk→∞

∫Ωukv =

∫Ωuv ∀v ∈ D(Ω)

Teorema A.2.2. Se uk → u nel senso delle distribuzioni e w e una distribuzione asupporto compatto allora uk ∗ w → u ∗ w nel senso delle distribuzioni (Cfr: [GA3, pag.340]).

Analisi Funzionale 89

Teorema A.2.3. Siano u e v ∈ L1loc((R)) e si supponga che sia soddisfatta una delle

seguenti due condizioni:i) v ha supporto compattoii) u, v hanno supporto incluso in [0,+∞)Allora per quasi ogni x reale la convoluzione e ben definita e appartiene a L1

loc((R)).

Dimostrazione. Per la dimostrazione si veda [GA3, pagg. 320-322].

Teorema A.2.4. Nelle ipotesi del teorema A.2.3 vale: supp(u ∗ v) = supp(u) + supp(v).

Dimostrazione. Per la dimostrazione si veda [GA3, pagg. 325-326].

Teorema A.2.5 (Cauchy-Lipschitz). Sia f : D → Rn, con D aperto di (R)n+1. Se:

1. f e continua in D.

2. f e localmente lipschitziana in D rispetto a y e uniformemente in t.

Allora per ogni punto (τ, ξ) ∈ D esiste un intorno Iδ di τ , Iδ = [τ − δ, τ + δ] nel quale edefinita una soluzione φ del problema di Cauchy (1.16) di [PS2]. Tale soluzione e unicanel senso che ogni altra soluzione coincide con φ nell’intervallo comune di definizione.

Dimostrazione. Per la dimostrazione si veda [PS2, pagg. 209-210]

A.2.1 Metodi di Linearizzazione per sistemi di EDO non lineari auto-nomi

L’analisi della stabilita degli equilibri di un sistema non lineare autonomo puo essereeffettuata sostituendo al sistema una sua approssimazione lineare nell’intorno dei puntidi equilibrio. Nonostante l’Analisi fornisca risultati globali per il sistema lineare, essendotale sistema una approssimazione di quello originale, le conclusioni che si traggono hannovalidita perlopiu locale. Dopo questa premessa introduciamo piu in dettaglio il metodo dilinearizzazione. Si consideri dunque un sistema autonomo del tipo:

y′ = f(y)

con f ∈ C1(Ω), Ω aperto di Rn e tale che f(yeq) = 0. L’ipotesi fondamentale per il seguitoe che la matrice Jacobiana di f in yeq, Df(yeq), sia non singolare. Essendo f ∈ C1(Ω) sipuo scrivere:

f(y) = Df(yeq)(y − yeq) + g(y − yeq)

con g(y − yeq) = o(||y − yeq||) per y → yeq. Poniamo A := Df(yeq) e consideriamo ilsistema:

z′ = Az

che prende il nome di sistema linearizzato in y = yeq. Il risultato fondamentale e chel’asintotica stabilita e/o l’instabilita del sistema linearizzato si puo trasferire al sistemaoriginale. Maggiori dettagli sono sintetizzati nel seguente risultato:

Teorema A.2.6. Si consideri un sistema differenziale omogeneo e autonomo del tipo:y′ = f(y). Allora:

90 Appendice

1. Se la matrice Jacobiana di f ha un autovalore con parte reale positiva in corrispon-tenza di un dato equilibrio y allora tale equilibrio e instabile.

2. Se la matrice Jacobiana di f ha tutti gli autovalori con parte reale negativa in corri-spontenza di un dato equilibrio y allora tale equilibrio e localmente asintoticamentestabile.

Dimostrazione. Per la dimostrazione si vedano i risultati sul metodo di linearizzazione in[PS2, pagg. 310-311]

Teorema A.2.7 (Poincare-Bendixon). Si consideri un sistema differenziale omogeneo eautonomo del tipo: y′ = f(y) definito su un insieme Ω ⊆ R2. Se un’ orbita Γ del sistemae interamente contenuta in un insieme limitato Ω′ ⊆ Ω allora si verifica soltanto una delleseguenti:

1. Γ e un punto di equilibrio

2. Γ tende a un punto di equilibrio

3. Γ e un’ orbita periodica

4. Γ tende ad un’ orbita periodica

(Cfr: [PS2, pag. 317]).

Elenco delle figure

1.1 Andamenti principali caratterizzanti il modello di Malthus . . . . . . . . . . . . . 31.2 Andamenti principali caratterizzanti il modello di Controllo Logistico . . . . . . . 41.3 Spazio delle fasi per il sistema di Lotka-Volterra classico . . . . . . . . . . . . . . 81.4 Soluzione del sistema di Lotka-Volterra classico in funzione del tempo . . . . . . . 91.5 Orbite di due sistemi di Lotka-Volterra classici chiarificanti il Teorema 1.4.3 . . . . 121.6 Soluzioni in funzione del tempo di due sistemi di Lotka-Volterra classici chiarificanti

il Teorema 1.4.3 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 121.7 Simulazione numerica atta a supportare la congettura 1.4.1 . . . . . . . . . . . . 151.8 Ingrandimenti della Figura 1.7 . . . . . . . . . . . . . . . . . . . . . . . . . . . 151.9 Spazio delle fasi per il sistema di Lotka-Volterra con memoria nel caso in cui

l’equilibrio (1.5.2b) sia attrattivo . . . . . . . . . . . . . . . . . . . . . . . . . . 181.10 Soluzioni in funzione del tempo del sistema di Lotka-Volterra con memoria nel caso

in cui l’equilibrio (1.5.2b) sia attrattivo . . . . . . . . . . . . . . . . . . . . . . 181.11 Spazio delle fasi per il sistema di Lotka-Volterra con memoria nel caso in cui

l’equilibrio (1.5.2c) sia attrattivo . . . . . . . . . . . . . . . . . . . . . . . . . . 191.12 Soluzioni in funzione del tempo del sistema di Lotka-Volterra con memoria nel caso

in cui l’equilibrio (1.5.2c) sia attrattivo . . . . . . . . . . . . . . . . . . . . . . . 19

2.1 Grafico del nucelo (2.2.1) nel caso base k = r = 1 . . . . . . . . . . . . . . . . . 252.2 Grafici del nucelo (2.2.1) al variare del PMN . . . . . . . . . . . . . . . . . . . . 262.3 Grafici del nucelo (2.2.1) al variare del PDN . . . . . . . . . . . . . . . . . . . . 262.4 Convoluzione tra il nucleo (2.2.1) e la funzione sin(t)H(t) quando il nucleo stesso

tenda alla delta di Dirac . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 312.5 Convoluzione tra il nucleo (2.2.1) e la funzione H(t) quando il nucleo stesso tenda

alla delta di Dirac . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 312.6 Convoluzione tra il nucleo (2.2.1) (r = k = 1) e la funzione sin(t)H(t) . . . . . . . 322.7 Convoluzione tra il nucleo (2.2.1) (r = 2, k = 1) e la funzione sin(t)H(t) . . . . . . 332.8 Convoluzione tra il nucleo (2.2.1) (r = 0.5, k = 1) e la funzione sin(t)H(t) . . . . . 332.9 Convoluzione tra il nucleo (2.2.1) (r = k = 1) e la funzione H(t) . . . . . . . . . . 342.10 Convoluzione tra il nucleo (2.2.1) (r = k = 1) e la funzione exp(t)H(t) . . . . . . 342.11 Convoluzione tra il nucleo (2.2.1) (r = k = 1) e la funzione exp(−t)H(t) . . . . . . 352.12 Convoluzione tra il nucleo (2.2.1) (r = k = 1) e la funzione t

10H(t) . . . . . . . . 35

3.1 Soluzione dell’equazione (3.1.1): repulsivita dell’equilibrio (3.2.1a) ed attrattivitadell’equilibrio (3.2.1b). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40

91

92 ELENCO DELLE FIGURE

3.2 Soluzione dell’equazione (3.1.1): repulsivita dell’equilibrio (3.2.1a) . . . . . . . . . 403.3 Soluzione dell’equazione (3.1.1): attrattivita dell’equilibrio (3.2.1b) . . . . . . . . 413.4 Soluzione dell’equazione (3.1.1): attrattivita dell’equilibrio (3.2.1b) . . . . . . . . 423.5 Soluzione dell’equazione (3.1.1): divergenza della soluzione per b < k . . . . . . . 433.6 Soluzione dell’equazione (3.1.1): divergenza della soluzione per b < k . . . . . . . 433.7 Soluzione dell’equazione (3.1.1): comportamento della soluzione per b = k . . . . . 443.8 Soluzione dell’equazione (3.1.1): considerazioni in relazione al Teorema 1.3.1 . . . 453.9 Soluzione dell’equazione (3.1.1): considerazioni in relazione al Teorema 1.3.1 ed al

Microteorema 3.3.1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 463.10 Soluzione dell’equazione (3.1.1): considerazioni in relazione al Teorema 1.3.1 ed al

Microteorema 3.3.1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 463.11 Soluzione dell’equazione (3.1.1): convergenza della soluzione a quella di un modello

di controllo logistico qualora r tenda a +∞ . . . . . . . . . . . . . . . . . . . . 483.12 Soluzione dell’equazione (3.1.1): convergenza della soluzione a quella di un modello

di controllo logistico qualora r tenda a +∞ . . . . . . . . . . . . . . . . . . . . 483.13 Andamenti degli autovalori della matrice Jacobiana associata all’equilibrio (3.2.1a) 503.14 Andamenti degli autovalori della matrice Jacobiana associata all’equilibrio (3.2.1b) 513.15 Soluzione dell’equazione (3.1.1): velocita di attrazione/repulsione degli equilibri al

variare del PDN . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 513.16 Soluzione dell’equazione (3.1.1): convergenza della soluzione a quella di un modello

di controllo logistico qualora il nucleo tenda alla massa di Dirac . . . . . . . . . . 523.17 Soluzione dell’equazione (3.1.1): convergenza della soluzione a quella di un modello

di controllo logistico qualora il nucleo tenda alla massa di Dirac . . . . . . . . . . 52

4.1 Soluzione del sistema (4.1.1): repulsivita dell’equilibrio (4.2.1a) e dell’equi-librio (4.2.1b) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58

4.2 Soluzione del sistema (4.1.1): repulsivita dell’equilibrio (4.2.1a) e dell’equi-librio (4.2.1b) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59

4.3 Soluzione del sistema (4.1.1): attrattivita dell’equilibrio (4.2.1b) . . . . . . . 604.4 Soluzione del sistema (4.1.1): attrattivita dell’equilibrio (4.2.1b) . . . . . . . 614.5 Soluzione del sistema (4.1.1): attrattivita dell’equilibrio (4.2.1b) . . . . . . . 614.6 Soluzione del sistema (4.1.1): attrattivita dell’equilibrio (4.2.1b) . . . . . . . 624.7 Autovalori della matrice Jacobiana associati all’equilibrio (4.2.1d): caso

sempre localmente asintoticamente stabile . . . . . . . . . . . . . . . . . . . 634.8 Autovalori della matrice Jacobiana associati all’equilibrio (4.2.1d): esistenza

di coppie r1, r2 che rendono l’equilibrio instabile . . . . . . . . . . . . . . . 644.9 Autovalori della matrice Jacobiana associati all’equilibrio (4.2.1d): regione

che rende l’equilibrio instabile . . . . . . . . . . . . . . . . . . . . . . . . . . 654.10 Soluzione del sistema (4.1.1): attrattivita dell’equilibrio (4.2.1d) nel caso in

cui non esista un coppia di valori r1, r2 che rende l’equilibrio instabile . . . 664.11 Soluzione del sistema (4.1.1): attrattivita dell’equilibrio (4.2.1d) nel caso in

cui esista un coppia di valori r1, r2 che renda l’equilibrio instabile . . . . . . 664.12 Soluzione del sistema (4.1.1): attrattivita dell’equilibrio (4.2.1d) nel caso in

cui non esistano dei valori del PDN che renda l’equilibrio instabile . . . . . 69

ELENCO DELLE FIGURE 93

4.13 Soluzione del sistema (4.1.1): attrattivita dell’equilibrio (4.2.1d) nel caso incui esista un valore del PDN che renda l’equilibrio instabile . . . . . . . . . 70

4.14 Soluzione del sistema (4.1.1): ciclo limite . . . . . . . . . . . . . . . . . . . . 744.15 Soluzione del sistema (4.1.1): ciclo limite . . . . . . . . . . . . . . . . . . . . 754.16 Soluzione del sistema (4.1.1): ciclo limite . . . . . . . . . . . . . . . . . . . . 754.17 Soluzione del sistema (4.1.1): ciclo limite . . . . . . . . . . . . . . . . . . . . 764.18 Soluzione del sistema semplificato dell’esempio in sezione 4.2.4: ciclo limite 764.19 Soluzione del sistema semplificato dell’esempio in sezione 4.2.4: ciclo limite 774.20 Soluzione del sistema (4.1.1): PDN tendenti a zero . . . . . . . . . . . . . . 804.21 Soluzione del sistema (4.1.1): PDN tendenti a zero . . . . . . . . . . . . . . 804.22 Soluzione del sistema (4.1.1): PDN tendenti a zero . . . . . . . . . . . . . . 814.23 Soluzione del sistema (4.1.1): PDN tendenti a zero . . . . . . . . . . . . . . 814.24 Soluzione del sistema (4.1.1): PDN tendenti a infinito . . . . . . . . . . . . 824.25 Soluzione del sistema (4.1.1): PDN tendenti a infinito . . . . . . . . . . . . 824.26 Soluzione del sistema (4.1.1): PDN tendenti a infinito . . . . . . . . . . . . 834.27 Soluzione del sistema (4.1.1): PDN tendenti a infinito . . . . . . . . . . . . 83

94 ELENCO DELLE FIGURE

Bibliografia

[BVIDE] T. A. Burton, Volterra Integral and Differential Equations, Mathematics inScience and Engineering, Volume 167, Academic Press, 1983.

[PS1] C.D.Pagani, S.Salsa, Analisi Matematica, Volume I, Masson, 1997.

[PS2] C.D.Pagani, S.Salsa, Analisi Matematica, Volume II, Masson, 1997.

[GA3] G.Gilardi, Analisi Tre, McGraw-Hill, 1994.

[NAL] W. K. Nicholson, Algebra Lineare, McGraw-Hill, 2002.

[BVNLCC] C.S. Bertuglia, F.Vaio, Non linearita, caos, complessita, Bollati Boringhieri,2003.

[DPRFA] F.Dercole, C.Piccardi, S.Rinaldi, Fondamenti di Automatica, retrieved fromhttp://www.elet.polimi.it, 2004.

[CACR] P. Colaneri, Appunti di Controllo Robusto, retrieved fromhttp://www.elet.polimi.it, 2006.

95

Nuclei di memoria in modelli biologici di Volterradigilander.libero.it/matteolesinigo/Matteo/Tesi di...

Documents

Transcript of Nuclei di memoria in modelli biologici di Volterradigilander.libero.it/matteolesinigo/Matteo/Tesi di...