Un'analisi sui sistemi di Reazione-Diffusione

Università degli Studi di Roma “La Sapienza”

Anno Accademico 2008-2009

Facoltà di Scienze Matematiche Fisiche e Naturali

Corso di Laurea Triennale in Matematica

Un’analisi sui sistemi di reazione-diffusione

Relatore Laureando

Prof. Andrea Terracina David Sarrocco

Matricola: 1153043 1

A Deborah,

Alla mia famiglia,

A Cinzia e Andrea,

Grazie..

Indice

Introduzione 3

1 Un Teorema di esistenza e unicità 6

1.1 Preliminari . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6

1.2 Teorema di esistenza e unicità . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12

2 Domini invarianti 16

2.1 Motivazioni . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16

2.2 Domini invarianti . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17

2.2.1 Una spiegazione “fisica” . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19

2.2.2 Il caso dei rettangoli . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20

3 Un teorema di confronto 24

3.1 Problemi associati . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24

3.2 Teorema di confronto . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29

4 Due esempi fisicamente interessanti 32

4.1 Modello Preda-Predatore . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32

4.1.1 Alla ricerca di un dominio invariante . . . . . . . . . . . . . . . . . . . . . . . . 33

4.1.2 Applichiamo il confronto . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35

4.2 Un Modello di Competizione

(scoiattoli rossi vs scoiattoli grigi) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36

4.2.1 Alla ricerca di un dominio invariante 2 . . . . . . . . . . . . . . . . . . . . . . . . 36

4.2.2 Applichiamo il confronto 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37

Bibliografia 39

2

Introduzione

Lo scopo che ci poniamo in questa tesi è quello di formulare un modello matematico, e una metodologia che permetta la descrizione dell’evoluzione di alcuni tra i più interessanti fenomeni che si presentano in natura. L’equazione che ci accompagna è detta di reazione-diffusione ed è della forma

∂u/∂t = D∆u + f(u), x ∈ Ω, t > 0 (1)

dove u ∈ R^m, D è una matrice m × m definita positiva (che verrà considerata diagonale), e f è una funzione regolare. Il laplaciano che compare nel secondo membro ci dà un’idea di come l’evoluzione dipenda dalla variazione della quantità u nello spazio; per tale motivo prende il nome di termine di “diffusione” in combinazione con la matrice D detta appunto matrice di diffusione.

Il termine f invece prende il nome di termine di “reazione” o termine “cinetico”, in quanto in assenza di diffusione la (1) diventa un’equazione differenziale ordinaria e la f ne definisce la cinetica.

L’importanza della (1) (e dei problemi di Cauchy associati) discende dal fatto che a seconda della matrice D e del termine di reazione si possono ottenere equazioni che descrivono le situazioni più differenti, dalla competizione tra specie animali alle applicazioni in finanza matematica, dai modelli preda-predatore alla dinamica degli enzimi. Per prima cosa cerchiamo di capire da dove esca fuori un’equazione di tale specie, ovvero preoccupiamoci di derivarne la forma.

Derivazione

Consideriamo u essere la densità di una determinata quantità misurabile (per esempio la densità di una popolazione in una determinata regione, o il livello di concentrazione di enzimi o proteine nel sangue), e supponiamo che la variazione di tale quantità/densità sia dovuta sia a un termine di “produzione locale” (per esempio la nascita/morte di elementi della popolazione), sia a un termine di trasporto/flusso (per esempio il trasferirsi da una regione all’altra della specie animale). Ovvero considerata una regione dello spazio Ω assumiamo valga la legge di bilancio:

d/dt ∫_Ω u dV = ∫_Ω f dV − ∫_∂Ω J · n dS, (2)

dove con f si intende appunto il termine di produzione/distruzione della quantità u, mentre J rappresenta il flusso uscente di tale quantità attraverso la superficie ∂Ω lungo la normale n.

3

INDICE 4

1

A questo punto se consideriamo Ω una regione abbastanza regolare e il flusso J di classe C¹(Ω) possiamo applicare il Teorema della divergenza

∫_∂Ω J · n dS = ∫_Ω div(J)dV (3)

e quindi sostituendo in (4) e supponendo l’uguaglianza ottenuta valga per ogni regione Ω si ottiene passando alla forma localizzata

∂u/∂t = f − div(J) (4)

Ora per quanto riguarda f ci occuperemo di funzioni dipendenti solo dallo stato u e non dal tempo (direttamente), mentre riguardo il flusso J si può considerare essere direttamente proporzionale alla variazione (componente per componente) della quantità u ovvero del suo gradiente. Quindi in formule

J = −d∇u (5)

dove d è una quantità che può dipendere da vari fattori, ma che noi considereremo costante e positiva. Infatti se la variazione è negativa (ovvero ho una diminuzione della quantità u) allora il flusso uscente dovrà necessariamente essere positivo (la quantità sta uscendo), se invece è positiva allora il flusso uscente sarà negativo, ovvero ho un’entrata di tale quantità in Ω. La (5) va sotto il nome di legge di Fourier o anche legge di Fick.

In finale quindi sostituendo l’espressione della (5) nella (4) si ottiene la forma scalare della (1) ovvero:

∂u/∂t = d∆u + f(u) (6)

Se la quantità u è vettoriale (u = u₁, ..., u_m) le varie componenti possono indicare ciascuna ad esempio un agente chimico presente in una regione dello spazio, o una specie animale presente su un territorio, in tal caso le formule con integrali vanno intese componente per componente e ciò che si sta considerando è il sistema

∂u₁/∂t = d₁∆u₁ + f₁(u₁, .., u_m)

..

∂u_m/∂t = d_m∆u_m + f_m(u₁, .., u_m) (7)

Ora, tornando all’obiettivo principe di questa trattazione (l’analisi di equazioni della forma (1)), facciamo un brevissimo quadro dei passi che verranno seguiti.

INDICE 5

In primo luogo dimostreremo che sotto opportune ipotesi l’equazione (1) associata a un particolare problema di Cauchy, ammette una soluzione unica perlomeno per tempi brevi (ovvero localmente nel tempo).

Fatto ciò però osserveremo che, dato il problema di Cauchy, qualora si fosse in grado di definire una limitazione a priori di quella che ne sarà la soluzione, allora si otterrà che tale soluzione esisterà globalmente nel tempo.

Introdurremo quindi il concetto di dominio invariante per un sistema di equazioni e mostreremo che se tale sistema ammette un dominio invariante allora la limitazione a priori è garantita.

Daremo quindi condizioni sufficienti (ma anche necessarie) affinché un’equazione-sistema ammetta tale dominio.

La nostra analisi poi ci porterà a utilizzare tale regioni per dimostrare un importante teorema che lega la soluzione della nostra equazione con la soluzione di un’equazione differenziale ordinaria, e ciò ci permetterà di ottenere informazioni per tempi grandi (quindi asintoticamente in t) sul comportamento dell’equazione o del sistema di equazioni.

Infine nell’ultima parte ci occuperemo di presentare due modelli classici nei quali tale analisi risulta esauriente e particolarmente interessante: il modello di competizione tra due specie animale, e il modello preda-predatore.

Capitolo 1

Un teorema di esistenza e unicità

1.1 Preliminari

L’obiettivo principale di questo capitolo è arrivare alla formulazione di un Teorema di Esistenza e Unicità per un problema di Cauchy per un particolare sistema di equazioni a derivate parziali di tipo parabolico. Le equazioni di cui ci occupiamo sono appunto le equazioni di reazione-diffusione.

In realtà la situazione che ci poniamo di studiare è abbastanza specifica, ovvero consideriamo il problema di Cauchy di un sistema in un dominio spaziale Ω ⊂ Rⁿ limitato con condizioni di Neumann al bordo. Ovvero siamo interessati a risolvere

∂u/∂t = c²∆u + f(u), (x, t) ∈ Ω × (t > 0)

u(x, 0) = u⁰(x)

∂u/∂n = 0, x ∈ ∂Ω × (t > 0) (1.1)

con l’uso della seguente notazione: c² = (c₁², c₂², .., c_m²)

u = ((u₁(x, t), u₂(x, t), ..., u_m(x, t))),

u⁰ = (u₁⁰(x, t), u₂⁰(x, t), ..., u_m⁰(x, t)),

f: R^m → R^m continua con derivata continua e tale che f(0) = 0,

u⁰ continuo con derivata continua e tale che soddisfi le condizioni al bordo,

con ‖u‖_∞ = max(‖u₁‖_∞, ..., ‖u_m‖_∞), ‖u‖_∞ = max_i max_x∈Ω (u_i(·, t)) (t > 0),

Ovvero stiamo considerando un sistema di m equazioni legate tra loro per il tramite del termine di reazione f.

In realtà diversi risultati che otterremo necessitano della sola lipschitzianità per la f (sia nel caso di sistemi, che nel caso scalare), per questo in alcuni punti della trattazione la considereremo solamente lipschitziana.

Per arrivare al nostro obiettivo procederemo per gradi, ovvero in primo luogo cercheremo la soluzione nel caso scalare e unidimensionale (x ∈ R considerando l’equazione omogenea (senza il termine di reazione),

6

1.1 Preliminari 7

termine di reazione), generalizzando poi quindi al caso di sistemi, infine tramite il metodo Duhamel definiremo la condizione che la soluzione del problema (1.1) deve verificare per essere tale. A quel punto potremo enunciare e dimostrare il teorema che pone le ipotesi per le quali tale condizione è verificata. La scelta di lavorare nel caso unidimensionale è per semplicità, i risultati ottenuti si estendono senza troppe difficoltà a più dimensioni.

Nel seguito indicheremo le derivate parziali con pedici sulla funzione.

Passo 1

Considerata l’equazione scalare e il problema associato

u_t = c²u_xx, x ∈ (0, π), (t > 0)

u(x, 0) = u⁰(x)

u_x(t) = 0, x ∈ {0, π} (t > 0) (1.2)

estendiamo il dato iniziale u⁰(x) per parità tra [−π, π] e considerata la base ortonormale di tale intervallo {1/√(2π), cos kx/√π, sin kx/√π, k = 1,2..}, supponendo che tale dato iniziale prolungato sia sviluppabile in serie di Fourier nella forma

u⁰(x) = 1/√(2π) ∫_−π^π u⁰(y) dy + Σ_k=1^∞ 1/√π ∫_−π^π u⁰(y)cos ky dy

si ottiene facilmente un’espressione per il dato iniziale (ricordando che ora u⁰(y) è pari) lo sviluppo della forma

u⁰(x) = 1/π [u₀⁰⁰ / 2 + Σ_k=1^∞ u_k⁰ cos kx]

con

u₀⁰⁰ = ∫₀^π u⁰(y)dy, u_k⁰ = ∫₀^π u⁰(y)cos(ky)dy

Lo sviluppo del dato iniziale ci suggerisce di cercare una “soluzione” della forma

u(x, t) = 1/π [a₀(t)/2 + Σ_k=1^∞ a_k(t)cos kx] (1.3)

A tale espressione, comunque, si giunge se si tenta di risolvere il sistema (1.2) per separazione di variabili e prendendo la serie della famiglia di soluzioni trovate. Ora considerando la serie (1.3) termine a termine si ottiene per ognuno di essi (sostituendolo nell’equazione, e ricordando che deve valere u(x, 0) = u⁰(x)) il problema ai coefficienti

ȧ_k(t) = −c²k²a_k(t), k = 0, 1, 2....

a_k(0) = u_k⁰ (1.4)

la cui soluzione è data da a_k(t) = u_k⁰e^−c2k2t se k = 0, 1, 2, ...

Allora sostituendo i parametri a_k(t) si ottiene una “soluzione” della forma

u(x, t) = 1/π [u⁰⁰₀/2 + Σ_k=1^∞ (∫₀^πu⁰(y)cos(ky)cos(kx)e^−c2k2tdy)] (1.5)

Ora, per costruzione ogni termine della serie soddisfa il problema di Cauchy (1.2), per cui ogni somma finita di tali termini lo soddisfa; verifichiamo che anche la serie lo soddisfa.

• Il dato iniziale è soddisfatto (per come è costruita la soluzione).
• Se potessimo derivare termine a termine sotto il segno di derivata allora, poiché ogni termine soddisfa il bordo e l’equazione, avremmo concluso.

Sia u_k il termine generico della serie, allora |u_k| ≤ ‖u⁰‖_∞e^−c2k2τ e la cui serie converge, per cui la nostra serie converge (totalmente) in [τ, ∞) ∀τ > 0.

A questo punto ci basta verificare che la serie delle derivate seconde rispetto a x e quelle delle derivate prime rispetto a t converge uniformemente:

|(u_xx(x, t))| = |−(2/π)Σ_k=1^∞(∫₀^πu⁰(y)cos(ky)cos(kx)c²k²e^−c2k2tdy)| ≤ ‖u⁰‖_∞c²k²e^−c2k2τ

la cui serie converge in [τ, T] con τ > 0, inoltre poiché tale stima è vera ∀t ∈ [τ, T] e ∀x ∈ (0, π) la serie converge totalmente e quindi uniformemente.

Analogamente vale ∀x ∈ (0, π) e ∀t ∈ [τ, T]

|(u_t(x, t))| = |−(2/π)Σ_k=1^∞(∫₀^πu⁰(y)cos(ky)cos(kx)k²e^−c2k2tdy)| ≤ ‖u⁰‖_∞k²e^−c2k2τ

Per cui è possibile effettuare la derivazione termine a termine e quindi l’equazione è verificata.

Con la stessa argomentazione si verifica che si può derivare una volta rispetto a x per cui è soddisfatta la condizione di Neumann.

Osservazione 1

Poiché la serie di u(x, t) converge uniformemente e poiché i termini u_n sono continui si ha che u(x, t) è continua. Inoltre poiché anche le serie derivate hanno la stessa proprietà, la soluzione risulta essere C^∞(0, π) per ogni t > 0, e continua nel tempo (per t > 0).

Verifichiamo infine la continuità nello spazio e nel tempo per t = 0, per farlo ci basta verificare la convergenza totale in [0,T] (quindi includendo lo 0) della serie (1.3). Vale infatti:

|a_k(t)cos kx| = |a_k(t)| = |a_k⁰e^−c2k2t| ≤ |a_k⁰| ≤ ‖u⁰‖_∞

1.1 Preliminari 9

1

e poiché stiamo assumendo che il dato iniziale sia su un insieme limitato, si verifica facilmente (con un’integrazione per parti) che la serie dei suoi coefficienti di Fourier converge totalmente in [0,T], ∀T > 0, per la nostra serie.

Osservazione importante

Se definiamo la funzione

g(x, y, t) = 1/π + (2/π)Σ_k=1^∞cos(kx)cos(ky)e^−c2k2t (1.6)

allora la soluzione (1.3) del nostro problema omogeneo (1.2) si può scrivere nella forma

u(x, t) = ∫₀^πg(x, y, t)u⁰(y)dy (1.7)

Occupiamoci di verificare due importanti proprietà della funzione g(x, y, t) definita in (1.6).

Proprietà 1

g(x, y, t) ∈ C^∞((0, π) × (0, π)) ∀t > 0

Dimostrazione. sup_x|e^−c2k2tcos kx cos ky| = e^−c2k2t che come già detto è il termine di una serie convergente per t > 0, da cui quindi già si osserva che g(x, y, t) è continua nelle variabili spaziali.

Inoltre considerata la serie delle derivate (rispetto a x e a y) si ottiene nuovamente lo stesso risultato, con la stessa maggiorazione, per cui la funzione è derivabile con derivata continua.

Iterando il processo si ha che g(x, y, t) ∈ C^∞((0, π) × (0, π)) per ogni tempo fissato.

Proprietà 2

∫₀^πg(x, y, t) = 1 ∀t > 0 fissato, (1.8)

dove il primo termine della (1.8) è da intendersi sia come integrale sulle y avendo fissato x, sia come integrale lungo le x avendo fissato y.

Dimostrazione.

∫₀^π∫₀^πg(x, y, t) = ∫₀^π∫₀^π[1/π + (2/π)Σ_k=1^∞e^−c2k2tcos(kx)cos(ky)] = 1,

in quanto i termini della sommatoria sono tutti nulli ∀k > 0, dal momento che sia integrando lungo le x sia lungo le y l’integrale del coseno tra [0, π] risulta nullo.

1.1 Preliminari 10

A questo punto se consideriamo il sistema di equazioni:

u_t = c²u_xx, x ∈ (0, π)

u(x, 0) = u⁰(x)

u_x(x) = 0, x ∈ {0, π}

con c² = (c₁², c₂², ...., c_m²), e u⁰(x) = (u₁⁰, ...u_m⁰), la soluzione sarà data da:

u(x, t) = ∫₀^πG(x, y, t)u⁰(y)dy

con G(x, y, t) = diag(g₁(x, y, t), ..., g_m(x, y, t)), dove g_i(x, y, t) = 1/π + (2/π)Σ_k=1^∞e^−c_i2k2tcos(kx)cos(ky).

Passo 2

Il secondo passo è quello che ci porta a definire quale deve essere la forma della soluzione del problema (1.1), ovvero quale equazione (integrale) deve soddisfare l’ipotetica soluzione qualora esista. Data la complessità di alcuni aspetti tecnici, riportiamo qui i risultati che ci interessano, rimandando ad altri testi per i dettagli.

Anche in questo caso continuiamo a considerare x ∈ R per semplicità, ripetendo nuovamente che le stesse considerazioni si possono estendere per x ∈ Rⁿ.

Con le notazioni introdotte consideriamo il sistema

u_t = c²u_xx + h(x, t), x ∈ (0, π), t > 0

u(x, 0) = u⁰(x)

u_x(x) = 0, x ∈ {0, π} (1.9)

La soluzione per tale problema è data dalla somma di due soluzioni rispettivamente dei problemi

u_t = c²u_xx, x ∈ (0, π), t > 0

u(x, 0) = u⁰(x)

u_x(x) = 0, x ∈ {0, π}