Metodi Matematici e

Statistici II

Versione Provvisoria

29 marzo 2001

La matematica e stata applicata ai problemi

reali fin dall’antichita

Il termine modello matematico e abbastanza

recente.

Si intende il tentativo di descrivere

un fenomeno mediante equazioni di-

sequazioni o altri strumenti matema-

tici, in modo da poter prevedere l’e-

voluzione del fenomeno stesso in di-

verse condizioni oppure in modo da

determinare una situazione ottimale.

Il processo di modellizzazione passa attraver-

so vari stadi

29 marzo 2001 1

Identificazione dei meccanismi e del-

le leggi che regolano il problema da

modellizzare: Formulare il problema

a parole e documentare i dati.

Formulazione delle ipotesi a parti-

re dalle quali il modello va costrui-

to: Individuare gli elementi importan-

ti e quelli che invece possono essere

ignorati.

Costruzione del modello attraverso la

traduzione della formulazione in lin-

guaggio matematico di quanto si e

gia definito in linguaggio naturale.

29 marzo 2001 2

Analisi del modello: Risolvere il pro-

blema matematico definito al punto

precedente.

Interpretazione del modello: Con-

frontare la soluzione del problema

matematico con la descrizione del

fenomeno che si vuole modellizzare.

Validazione del modello: confronto

tra le previsioni e i dati sperimentali.

Implementazione del modello: uso

del modello per prevedere l’evoluzio-

ne del fenomeno in tempi futuri o di-

verse condizioni. Determinazione di

parametri ottimali

29 marzo 2001 3

Modelli differenziali

di

crescita delle popolazioni

29 marzo 2001 4

Modelli di crescita per popolazioni

isolate,

comprendono

• crescita dell’inflazione,

• calcolo dell’interesse composto,

• decadimento di materiale radioattivo

(modello esponenziale modello logistico)

29 marzo 2001 5

Modelli di crescitaper popolazioni coesistenti

• in competizione

• in cooperazione.

• in lotta.

comprendono

• diffondersi di un’epidemia,

• combattimento tra due eserciti

29 marzo 2001 6

Modello di crescita esponenziale

E il piu semplice modello di crescita di una

popolazione

e stato proposto da Malthus nel 1798

prevede che la crescita di una popolazione sia

proporzionale al numero di individui.

Se il numero di individui di una

popolazione al tempo t e

n(t)

la variazione n(t) della popolazione e

assunta proporzionale al numero n(t)

di individui secondo una costante c

(tasso di crescita).

29 marzo 2001 7

Si avra n(t) = cn(t)

n(t0) = n0(1)

e si ottiene facilmente che

n(t) = n0ec(t−t0) (2)

Per determinare la costante c occorrono due

valori

n0 = n(t0) n1 = n(t1)

dalla 2 per t = t1

n1 = n0ec(t1−t0) e c =

ln(n1n0

)

t1 − t0(3)

29 marzo 2001 8

Nel caso in cui n1 = 2n0 si parla di tempo di

raddoppio mentre nel caso in cui n1 = 12n0 si

parla di tempo di dimezzamento.

Per come e stato definito, c rappresenta il

tasso di crescita istantaneo

c =n(t)

n(t)

In molti casi puo essere comodo considerare

il tasso di crescita mediato su un periodo T .

29 marzo 2001 9

Se c e il tasso di accrescimento della quantita

n0 nel periodo T (cioe e la frazione di n0

che si va ad aggiungere ad n0) avremo che il

valore finale e

nT = n0 + cn0 = n0(1 + c) (4)

In altre parole l’accrescimento dipende dal

valore iniziale n0 e dal tasso di crescita.

Possiamo usare anche un modello discreto.

Il tasso di accrescimento viene assunto co-

stante su un certo intervallo di tempo

29 marzo 2001 10

Si consideri il caso di una popolazione

di individui il cui tasso di accrescimen-

to medio su un periodo T e γ.

Se assumiamo l’accrescimento sull’in-

tero periodo solo gli individui presenti

all’istante iniziale

contribuiscono all’aumento della po-

polazione;

Possiamo pero supporre di misurare

l’accrescimento piu volte nel periodo

T a intervalli regolari di lunghezza Tk ,

Si verifica che la popolazione, do-

po un tempo T , sara cresciuta an-

che per opera degli individui nati nei

periodi [t0 + hTk , t0 + (h + 1)Tk ] per

h = 0,1, .., k − 1.29 marzo 2001 11

Si consideri un capitale iniziale n0 che

assicura un rendimento γ nel periodo

T .

Calcolare il rendimento sull’intero pe-

riodo T significa supporre di contabi-

lizzare gli interessi solo alla fine del-

l’intero periodo (interesse semplice)

Mentre frazionare il periodo signifi-

ca capitalizzare gli interessi su ognu-

no dei sottoperiodi di lunghezza Tk

e quindi percepire anche gli interessi

sugli interessi (interesse composto).

29 marzo 2001 12

In entrambi i casi si ottengono modelli discre-

ti che sono governati da regole di ricorren-

za che identificano una successione di valori

reali.

Se il tasso di crescita istantaneo e c, nel pe-

riodo T avremo un tasso di crescita pari a

γ = cT e la quantita in crescita (o decrescita

se c < 0) puo essere calcolata mediante la

n(t0 + T ) = n(t0) + cTn(t0)

n(t0) = n0(5)

Mentre sull’intero periodo, se consideriamo k

suddivisioni, si ha

n(t0) = n0

29 marzo 2001 13

e γk = cT

k e il tasso di accrescimento medio su

[hTk , (h+ 1)Tk ]

n

(t0 + h

T

k

)=

n

(t0 + (h− 1)

T

k

)+cT

kn

(t0 + (h− 1)

T

k

)(6)

per 1 ≤ h ≤ k

Si vede subito che

n

(t0 + h

T

k

)=(1 +

cT

k

)n

(t0 + (h− 1)

T

k

)=

=(1 +

cT

k

)2n

(t0 + (h− 2)

T

k

)=

=(1 +

cT

k

)hn(t0) =

=(1 +

cT

k

)hn0 (7)

29 marzo 2001 14

Confrontando le soluzioni della 2 e della 7,

avremo, dopo un tempo T

n(t0 + T ) = n0ecT o n(t0 + T ) = n0

(1 +

cT

k

)k(8)

Si rileva che per k → +∞

(1 +

cT

k

)k→ ecT

29 marzo 2001 15

Il modello di crescita Malthusiano prevede

un incremento della popolazione che non e

ragionevole.

Nel 1960 la popolazione mondiale era di 3×109 individui; assumendo un fattore di cre-

scita del 2% per anno si ha c =ln(

n1n0

)

t1−t0 =

ln1.02 ≈ .0198

n(t) = 3× 109eln(1.02)(t−1960)

per cui

1960 1970 19803× 109 3.7× 109 4.5× 109

1990 2000 20205.4× 109 6.6× 109 9.8× 109

2040 2060 208014.6× 109 21.7× 109 32.2× 109

2100 2200 230047.9× 109 347.4× 109 2516.5× 109

29 marzo 2001 16

Il modello esponenziale permette anche di

rappresentare, la rendita di un investimento

a tasso fisso o gli effetti dell’inflazione.

Ad esempio sapendo che il prezzo di un de-

terminato bene e passato da 1000 lire 20 an-

ni orsono a lire 10000 attualmente, possia-

mo calcolare il tasso di inflazione ricavando

c nell’equazione

10000 = 1000ec20

da cui

c =ln(10)

20= 0.115

per ottenere che nel periodo il tasso di infla-

zione e stato di circa l’11%.

29 marzo 2001 17

Anche il decadimento radioattivo e modelliz-

zato per mezzo di una legge esponenziale.

Ogni elemento radioattivo decade in maniera

proporzionale alla propria massa.

Per descrivere il fenomeno basta conoscere il

tempo necessario perche la massa si riduca a

meta (tempo di dimezzamento.)

Ad esempio il torio ha un tempo di dimezza-

mento di 1.65× 1010 anni;

pertanto

q = q0ect

dove c si puo ricavare mediante la

c =ln(1

2)

1.65× 1010= −1.33× 10−18

29 marzo 2001 18

Modello di crescita logistica

La crescita di una popolazione non e ben rap-

presentata da un semplice modello esponen-

ziale.

Si trascura completamente l’effetto dell’af-

follamento e della saturazione dell’ambiente.

Una popolazione in crescita consuma risorse

dall’ambiente e tali risorse non sono illimitate

e ne condizionano lo sviluppo.

Un modello che tenga conto anche di que-

ste condizioni fu proposto da Verhulst , ed

assume un fattore di accrescimento del tipo

a− bn a, b > 0

dove n e il numero di individui.

Il tasso di crescita diminuisce con il nume-

ro degli individui secondo una costante di

proporzionalita b

29 marzo 2001 19

Il tasso di crescita e‘ lineare decrescente

n =a

b

e il massimo numero di individui sostenibili

dall’ambiente.

Il numero di individui n e regolato dal seguen-

te problema di Cauchy:

n(t) = an(t)− bn2(t) = n(t)(a− bn(t))

n(t0) = n0

(9)

prende il nome di equazione logistica,

e possibile trovare informazioni sulla soluzio-

ne (funzione logistica) senza risolverla.

29 marzo 2001 20

Ci sono due punti di equilibrio (solu-

zioni costanti) in corrispondenza dei

valori

n(t) = 0 n(t) =a

b= nmax

Corrispondono al caso in cui la po-

polazione si e estinta o ha saturato

l’ambiente in cui vive.

Il modello esponenziale prevede una crescita

illimitata ed ammette come unico punto di

equilibrio l’estinzione della popolazione,

La crescita logistica prevede anche

un valore di equilibrio non nullo che

abbiamo indicato con nmax

29 marzo 2001 21

poiche an− bn2 > 0 per 0 < n < ab La

soluzione e crescente se n0 <ab ed e

decrescente se n0 >ab

inoltre, poiche

n(t) = (a− 2bn(t))n(t)

La soluzione e convessa ( n crescen-

te), se n(t) < a2b ed e concava ( n

decrescente), per ab > n(t) > a

2b.

per n > ab la soluzione risulta

convessa.

29 marzo 2001 22

Andamento della popolazione

nel caso di crescita logistica

a = 4 b = 1

29 marzo 2001 23

La soluzione tende ad un valore asin-

totico nmax,

crescendo se n0 < nmax

decrescendo se n0 > nmax

(rimane costante se n0 = nmax).

Per n(t) = a2b = nmax

2 c’e un punto di

flesso

La crescita e molto rapida fino a

che n(t) = a2b = nmax

2 rallenta

avvicinandosi all’asintoto.

29 marzo 2001 24

L’equazione logistica

n(t) = an(t)− bn2(t) = n(t)(a− bn(t))

n(t0) = n0

(10)

puo essere facilmente integrata separando le

variabili.

Ammette sempre una ed una sola soluzione

Ammette due soluzioni costanti

n(t) = 0 , n(t) =a

b(11)

(soluzioni di equilibrio del modello)

identificano le situazioni in cui il numero di

individui e stabile29 marzo 2001 25

Possiamo determinare la soluzione come se-gue

n(t)

n(t)(a− bn(t))= 1

∫ tt0

n(s)

n(s)(a− bn(s))ds = t− t0

∫ n(t)n(t0)

dx

x(a− bx)= t− t0 (12)

e

1

x(a− bx)=

1

a

(1

x+

b

a− bx

)da cui∫ n(t)

n0

dx

x(a− bx)dx =

1

a

(ln |x|−ln |a−bx|

)n(t)n0

=

=1

aln

(n(t)

a− bn(t)

)(a− bn0

n0

)(13)

29 marzo 2001 26

Ne viene

t− t0 =1

aln c

(n(t)

a− bn(t)

)

dove n(t)a−bn(t) ha lo stesso segno di c = n0

a−bn0in ciascuno degli intervalli in cui e definita

Si ricava

n(t) =a

b+ ce−a(t−t0)=

1

β + γe−a(t−t0)(14)

dove β = ba e γ = c

a e le costanti a, β e γ

possono essere determinate se si conoscono

n(t0 − h) = p0 (15)

n(t0) = p1 (16)

n(t0 + h) = p2 (17)

29 marzo 2001 27

infatti si ha

1

n(t)= β + γe−a(t−t0) (18)

ed e sufficiente risolvere il sistema

1p0

= β + γeah

1p1

= β + γ1p2

= β + γe−ah(19)

posto H = e−ah

1p0

= β + γ 1H

1p1

= β + γ1p2

= β + γH

(20)

29 marzo 2001 28

da cui

1p0− 1p1

= γ(

1H − 1

)= γ1−H

H1p1− 1p2

= γ(1−H)(21)

e

1p0− 1p1

1p1− 1p2

=1

H(22)

se ne deduce che

1

H=p1p2(p1 − p0)

p0p1(p2 − p1)H =

p0p1(p2 − p1)

p1p2(p1 − p0)(23)

e poi facile ricavare γ e β dalle 20,21.

29 marzo 2001 29

Modello di

crescita logistica

con prelievo costante

Supponiamo di esaminare una popolazionecon crescita logistica sottoposta ad un pre-lievo costante nel tempo.

Ad esempio un allevamento ittico, da cui vie-ne pescata una quantita fissa di pesce.

Oppure una specie in ambiente protetto esenza antagonisti il cui sviluppo deve esse-re controllato mediante un prelievo periodicoper mantenere il numero di esemplari stabile.

La popolazione si sviluppa con tasso di cre-scita a− bn

Subisce una costante diminuzione pari ad unaquantita c.

Il valore di c e critico e puo portare all’estin-zione della popolazione.

29 marzo 2001 30

E utile determinare il massimo prelievo che

non causa l’estinzione della popolazione

occorre modificare il modello logistico.

n(t) = an(t)− bn2(t)− c

n(t0) = n0(24)

In ogni caso avremo che a, b, c > 0

L’equazione e a variabili separabili e puo es-

sere integrata.

Ammette sempre una ed una sola soluzione

Ammette soluzioni costanti nel caso in cui

l’equazione algebrica

−bx2 + ax− c = 0

abbia soluzioni reali.

29 marzo 2001 31

Possiamo distinguere tre casi

a2 − 4bc < 0

Non ci sono soluzioni costanti e si ha

n(t)

an(t)− bn2(t)− c= 1 (25)∫ t

t0

n(s)

an(s)− bn2(s)− cds = t− t0 (26)∫ n(t)

n(t0)

dx

ax− bx2 − c= t− t0 (27)

29 marzo 2001 32

si ottiene

1

ax− bx2 − c= −

1

b

1(x− a

2b

)2− a2

4b2+ c

b

=

= −1

b

1(x− a

2b

)2+ δ

=

= −1

bδ

1(1√δ

(x− a

2b

))2+ 1

= (28)

dove si e posto

δ =4bc− a2

4b2

29 marzo 2001 33

da cui∫dx

ax− bx2 − c=

= −1

b√δarctan

1√δ

(x−

a

2b

)+ cost. (29)

e

−1

b√δarctan

1√δ

(n(t)−

a

2b

)+

+1

b√δarctan

1√δ

(n0 −

a

2b

)= t− t0 (30)

inoltre

n(t) =a

2b+

+√δ tan

(t0 +

1

b√δarctan

1√δ

(n0 −

a

2b

)− t

)(31)

29 marzo 2001 34

Andamento della popolazione

nel caso di crescita logistica

con prelievo costante nel caso a2 − 4bc < 0

a = 2 b = 1 c = 4 a2−4bc = −12 < 0

29 marzo 2001 35

a2 − 4bc = 0

Esiste una soluzione doppia per l’equazione

−bx2 + ax− c = 0

e quindi una soluzione costante

n(t) = α =a

2b

Si ha inoltre

1

ax− bx2 − c= −

1

b

1(x− a

2b

)2 (32)

da cui∫ n(t)n0

dx

ax− bx2 − c=

1

b

1(n(t)− a

2b

) −−

1

b

1(n0 − a

2b

) = t− t0 (33)

e

n(t) =1

2b

(a+

2bn0 − a

1 + (n0b− a2)(t− t0)

)(34)

29 marzo 2001 36

Andamento della popolazione

nel caso di crescita logistica

con prelievo costante nel caso a2 − 4bc = 0

a = 2 b = 1 c = 1 a2 − 4bc = 0

29 marzo 2001 37

a2 − 4bc > 0

−bx2 + ax− c = 0

ammette due soluzioni reali e distinte

α =a

2b+

√a2 − 4bc

4b2, β =

a

2b−

√a2 − 4bc

4b2

e ci sono due soluzioni costanti (di equilibrio

per il sistema).

Si ha

−bx2 + ax− c = −b(x− α)(x− β)

e

29 marzo 2001 38

1

−bx2 + ax− c=

b

α− β

(1

x− α−

1

x− β

)(35)

da cui∫ 1

−bx2 + ax− cdx =

=b

α− β(ln |x− α| − ln |x− β|) + ϕ(x) (36)

e ϕ e una funzione costante a tratti∫ n(t)n0

1

−bx2 + ax− cdx =

=b

α− β

(lnn(t)− α

n(t)− β

n0 − β

n0 − α

)= t− t0 (37)

da cui

n(t) =α− βn0−α

n0−βe(t−t0)α−βb

1− n0−αn0−βe

(t−t0)α−βb(38)

29 marzo 2001 39

Andamento della popolazione con prelievo costante

nei casi in cui a2 − 4bc > 0

a = 5 b = 1 c = 4 a2−4bc = 9 > 0

29 marzo 2001 40

Andamento della popolazione con prelievo costante

nei casi in cui a2 − 4bc > 0

a = 3 b = 1 c = 2 a2−4bc = 1 > 0

29 marzo 2001 41

Come si puo facilmente osservare dal graficodella soluzione n(t) il modello rende ragionedei seguenti fatti

1. se a2− 4bc < 0 la popolazione si estinguein un tempo finito

2. se a2 − 4bc = 0

• la popolazione si estingue in un tempofinito se n0 <

a2b

• mentre tende a stabilizzarsi sul valorea2b asintoticamente se n0 >

a2b

3. se a2 − 4bc > 0

• la popolazione si estingue in un tempofinito se n0 < β

• mentre tende a stabilizzarsi sul valoreβ asintoticamente se n0 > β

29 marzo 2001 42

Per c = a2

4b risulta a2 − 4bc = 0

c e la piu alta quantita che si puo prelevare

senza causare l’estinzione in un tempo finito

Tale quantita viene detta rendita massima.

In presenza di un prelievo costante il valore

asintotico della popolazione diminuisce

da ab ad α = a

2b +

√a2−4bc

4b2< a

b .

Per a2 − 4bc = 0 si ha che α = a2b

29 marzo 2001 43

Modelli di crescita

con prelievo non costante

Nel modello di crescita logistica il tasso di

crescita della popolazione x e dato da a− bx

Avremo che

x(t) = (a− bx(t))x(t) (39)

Se il prelievo e costante, il tasso di crescita

e a− bx− h

e l’equazione

x(t) = (a− bx(t))x(t)− h (40)

29 marzo 2001 44

Se il prelievo e proporzionale al numero di

individui presenti. il termine che indica il

prelievo sara del tipo αx

e l’equazione diventa

x(t) = (a− bx(t))x(t)− αx(t) =

= ((a− α)− bx(t))x(t) (41)

non ce differenza rispetto al modello di cre-

scita logistica standard cambia solo il tasso

da a in a− α.

Si possono considerare tassi di crescita piu

generali

29 marzo 2001 45

ad esempio

f(x) =

(τ +

(β − x)x

1 + αx

)(42)

(f si annulla in x = 0 ed x = β ha un massimo

in√

1+αβ−1α ∈ (0, β))

che corrisponde ad una variazione della po-

polazione

x

(τ +

(β − x)x

1 + αx

)(43)

o piu semplicemente si possono considerare

a e b funzioni di t, da cui

x(t) = (a(t)− b(t)x(t))x(t) (44)

Equazione di Bernoulli, puo essere integrata

usando la sostituzione y(t) = 1x(t)

29 marzo 2001 46

Un modello di gestione

delle risorse ittiche

La crescita di una popolazione di pesci, puo

essere descritta mediante un modello logisti-

co del tipo

x(t) = f(x(t))− νy(t)x(t) (45)

dove f(x) = ax− bx2

e dove il termine νy(t)x(t) misura il prelievo

dovuto alla pesca proporzionale

al numero di individui x(t),

e alle risorse impegnate nella pesca y(t)

29 marzo 2001 47

Per avere una migliore rispondenza al mo-

dello reale possiamo considerare la crescita

regolata da una funzione

f(x) = x

(τ +

(β − x)x

1 + αx

)

Tuttavia supporremo solo che f sia positiva

e che

d

dx

(f(x)

x

)R 0 per x Q x (46)

Osserviamo che

d

dx

(f(x)

x

)=xf ′(x)− f(x)

x2

per cui

xf ′(x)− f(x) R 0 ⇔ x Q x

29 marzo 2001 48

Il ricavato della pesca e proporzionale al pe-

scato.

Si puo stimare nella misura y(t)(νpx(t)− c)

p e proporzionale al prezzo unitario di vendita

cy(t) e il costo delle risorse impiegate.

L’incremento delle risorse impegnate e pro-

porzionale a αy(t)(νpx(t)− c)

per cui

y(t) = αy(t)(νpx(t)− c) (47)

Quindi

x(t) = f(x(t))− νy(t)(t)x(t)

y(t) = αy(t)(νpx(t)− c)(48)

29 marzo 2001 49

E interessante studiare la stabilita del siste-

ma.

Il sistema ha un punto di equilibrio in

x =c

νpy =

f(x)

νx(49)

Il sistema linearizzato e

(x(t)− x)′ =f ′(x)(x(t)− x)− νy(x(t)− x)−

− νx(y(t)− y)(y(t)− y)′ = + αyνp(x(t)− x) +

α(y(t)− y)(νpx− c)(50)

29 marzo 2001 50

Posto ξ = x− x , η = y− y, il sistema diventa

ξ(t) = f ′(x)ξ(t)− νyξ(t)− νxη(t)

η(t) = αη(t)(νpx− c) + αyνpξ(t)(51)

la cui matrice dei coefficienti e(f ′(x)− νy −νxανpy 0

) f ′(x)− f(x)x −c

p

αpf(x)x 0

(52)

29 marzo 2001 51

Si ha

detA = αcνy > 0 (53)

quindi se la traccia di A che e f ′(x) − f(x)x e

positiva, allora la matrice ha autovalori con

parte reale positiva ed il sistema non e stabile

mentre se la traccia di A che e f ′(x)− f(x)x e

negativa, allora la matrice ha autovalori con

parte reale negativa ed il sistema e stabile

Se la traccia di A che e f ′(x) − f(x)x e nulla,

allora la matrice ha autovalori nulli; si vede

facilmente che la soluzione e stabile per il

sistema linearizzato; non si puo concludere

per il sistema originale.

29 marzo 2001 52

Per le ipotesi fatte f ′(x)− f(x)x R 0 se e solo

se x Q x;

infatti si ha f ′(x)− f(x)x e per la 46 f ′(x)− f(x)

xe positiva prima e negativa dopo x

Nel caso in cui il sistema sia instabile si puo

provare che per x − x piccolo esiste un ciclo

limite.

29 marzo 2001 53

Un modello di

ottimizzazione del profitto

Consideriamo una societa il cui scopo e quel-

lo di assicurare il massimo utile ai soci.

Supponiamo che gli utili siano proporzionali

al capitale investito

Se ogni utile e distribuito ai soci non e possi-

bile ottenere un incremento dell’utile stesso.

Se si reinveste parte degli utili per aumentare

il capitale si puo invece far aumentare l’utile

successivo.

29 marzo 2001 54

Il problema e

Determinare quale parte del profitto

deve essere distribuita subito in modo

che gli utili distribuiti ai soci su un

periodo fissato siano massimi

Il problema e di natura discreta, tuttavia e

comodo supporre che capitale ed utili siano

una funzione continua e derivabile del tempo.

Sia x(t) il capitale di cui dispone la societa

al tempo t;

Se il rendimento istantaneo e a in [t, t+ ∆t]

gli utili ammonteranno ad ax(t)∆t.

29 marzo 2001 55

Sia k ∈ [0,1], esprimiamo l’incremento di ca-

pitale ∆x nel periodo [t, t + ∆t] mediante

la

∆x = kax∆t (54)

Similmente il dividendo distribuito y(t) e

∆y = (1− k)ax∆t (55)

Da cui

x(t) = kax(t)

y(t) = (1− k)ax(t)

x(0) = α

y(0) = 0

(56)

dove α e il capitale iniziale

29 marzo 2001 56

Si tratta di un modello che richiama il mo-

dello di crescita malthusiana (esponenziale).

Si ottiene quindi che

x(t) =αeakt (57)

y(t) =

(1−k)α

k (eakt − 1) 0 < k ≤ 1

aαt k = 0(58)

Se k = 0 si ha redistribuzione dell’intero utile

il capitale rimane costante ed uguale ad α

Se k = 1 non si ha distribuzione di utile ed il

capitale cresce con legge esponenziale.

29 marzo 2001 57

Assegnato un periodo di tempo fissa-

to T come deve essere scelto k in mo-

do che sia massimo l’utile distribuito

al tempo T ?

Bisogna determinare per quali valori di k si

ottiene il massimo di

y(T ) =(1− k)α

k(eakT − 1) con 0 < k < 1

(59)

Posto ξ = akT la funzione da minimizzare e

φ(ξ) = αaT − ξ

ξ(eξ − 1) con 0 < ξ < aT

(60)

29 marzo 2001 58

Avremo

limξ→0

φ(ξ) = αaT limξ→aT

φ(ξ) = 0 (61)

inoltre possiamo calcolare che

φ′(ξ) =aTeξ

ξ2

(−1 + ξ −

ξ2

aT+ e−ξ

)(62)

occorre conoscere il segno di φ′ e l’eventuale

punto in cui si annulla.

E sufficiente studiare

ψ(ξ) =

(−1 + ξ −

ξ2

aT+ e−ξ

)(63)

29 marzo 2001 59

ed e utile calcolare

ψ′(ξ) = 1−2ξ

aT− e−ξ (64)

ψ′′(ξ) = −2

aT+ e−ξ (65)

Osserviamo ancora che

limξ→0

ψ′(ξ) = 0 limξ→aT

ψ′(ξ) = −1− e−aT < 0

(66)

limξ→0

ψ(ξ) = 0 limξ→aT

ψ(ξ) = −1 + e−aT < 0

(67)

29 marzo 2001 60

Possiamo distinguere due casi

2aT ≥ 1

dalla 65

ψ′′ ≤ 0

ψ′ e decrescente

per la 66 ψ′ ≤ 0

per la 64 e la 67 anche ψ ≤ 0

φ′ ≤ 0

φ sempre decrescente.

maxξ≥0

φ(ξ) = φ(0)

da cui k = 029 marzo 2001 61

2aT < 1

La 65 si annulla in un punto ξ0 = ln aT2 e

positiva prima e negativa dopo tale punto.

ψ′ ha un massimo in ξ0 e si annulla (per le

66) in un punto ξ1 < aT .

ψ e crescente prima di ξ1 e decrescente dopo,

per cui ψ ammette un punto di massimo in

ξ1 e si annulla (per le 67) in un punto ξ2,

ξ0 < ξ1 < ξ2 < aT .

φ ammette un punto di massimo in

ξ2 ∈ (lnaT

2, atT )

29 marzo 2001 62

29 marzo 2001 63

29 marzo 2001 64

Modelli di crescita

di due specie

in competizione

Consideriamo due specie x e y che condivi-

dono le risorse di uno stesso territorio

E ragionevole supporre che

Ciascuna specie abbia un tasso di crescita di

tipo logistico nel caso viva in un ambiente

isolato, cioe

x(t) = (a−Ax(t))x(t)

y(t) = (b−By(t))y(t)

x(t0) = x0y(t0) = y0

(68)

29 marzo 2001 65

Se si tiene conto dei mutui effetti

il tasso di crescita; della popolazione x e

a−Ax− αy

(diminuisce proporzionalmente a y)

il tasso di crescita della popolazione y e

b−By − βx.

(diminuisce proporzionalmente a x)

Ciascuna specie si sviluppa sottraendo risorse

all’altra.

29 marzo 2001 66

Il sistema che descrive lo sviluppo delle due

popolazioni diventa

x(t) = (a−Ax(t)− αy(t))x(t)

y(t) = (b−By(t)− βx(t))y(t)(69)

con

x(t0) = x0y(t0) = y0

(70)

Ove a, b, A,B, α, β sono positivi.

29 marzo 2001 67

Dalla prima equazione, dividendo per x ed

integrando otteniamo che

x(t)

x(t)= (a−Ax(t)− αy(t))

∫ x(t)x0

ds

s=∫ tt0(a−Ax(s)− αy(s))ds

lnx(t)

x0=∫ tt0(a−Ax(s)− αy(s))ds

e

x(t) = x0e

∫ tt0

(a−Ax(s)−αy(s))ds(71)

29 marzo 2001 68

In maniera simile otteniamo che

y(t) = y0e

∫ tt0

(b−By(s)−βx(s))ds(72)

Se ne deduce che

x(t), y(t) > 0

Dal sistema 69, per k = max{a, b} che

(x+ y)′ = (a−Ax(t)− αy(t))x(t) +

+ (b−By(t)− βx(t))y(t) ≤ k(x+ y) (73)

29 marzo 2001 69

Quindi x ed y sono limitate, in quanto

0 ≤ (x+ y) ≤ ekt ≤ ekT

e quindi esistono in ogni intervallo [0, T ].

Il sistema 69 ammette soluzioni costanti che

si ottengono risolvendo il sistema algebrico

0 = (a−Ax− αy)x

0 = (b−By − βx)y(74)

Le soluzioni costanti sono individuate dai pun-

ti di intersezione degli assi con le rette

R1 : (a−Ax− αy) = 0 (75)

R2 : (b−By − βx) = 0 (76)

29 marzo 2001 70

La retta R1 interseca gli assi nei punti

(ξ,0) , (0, η)

dove

ξ =a

Aη =

a

α

La retta R2 interseca gli assi nei punti

(ξ,0) , (0, η)

dove

ξ =b

βη =

b

B

Le rette R1 ed R2 si intersecano nel punto

(x, y)

dove

x =aB − αb

AB − αβ, y =

Ab− aβ

AB − αβ

29 marzo 2001 71

Pertanto le soluzioni costanti sono individua-

te dai punti

(0,0) , (0, η) , (ξ,0) , (x, y)

(77)

L’ultimo punto si considera solo nel caso in

cui

1. R1 ed R2 non siano parallele

(AB − αβ 6= 0)

2. si trovi nel primo quadrante

(x > 0 , y > 0).

Possiamo disegnare il campo di vettori asso-

ciato al sistema e studiare le orbite.

29 marzo 2001 72

Si presentano quattro casi nei quali

1. Prevale la popolazione x

2. Prevale la popolazione y

3. A seconda dei valori iniziali prevale x op-

pure y

4. Comunque si scelgano i valori iniziali, la

situazione tende ad un unico punto di

equilibrio

Nel primo, secondo e quarto caso c’e un pun-

to asintoticamente stabile

Nel terzo ci sono due punti asintoticamente

stabili ed un punto instabile

29 marzo 2001 73

Andamento del numero di individui

di due popolazioni in competizione

nel caso in cui la popolazione prevalente

dipenda dai dati iniziali

29 marzo 2001 74

Andamento del numero di individui

di due popolazioni in competizione

nel caso in cui

la popolazione prevalente sia la seconda

29 marzo 2001 75

Andamento del numero di individui

di due popolazioni in competizione

nel caso in cui le due popolazioni

tendano a raggiungere un equilibrio stabile

29 marzo 2001 76

Andamento del numero di individui

di due popolazioni in competizione

nel caso in cui la popolazione prevalente sia la prima

29 marzo 2001 77

Per studiare la stabilita delle soluzioni costan-

ti possiamo considerare il sistema linearizza-

to relativamente ad ognuna delle soluzioni

costanti (x, y).

Il sistema linearizzato in (x, y) e associato alla

matrice

M =

(a− 2Ax− αy −αx

−βy b− 2By − βx

)(78)

Per (x, y) = (0,0)

M =

(a 00 b

)(79)

M ha autovalori positivi e quindi la relativa

soluzione e instabile

29 marzo 2001 78

Per (x, y) = (ξ,0)

M =

(a− 2A a

A −α aA0 b− β aA

)=

=

(−a −α aA0 β

(bβ −

aA

)) =

=

(−a −α aA0 β(ξ − ξ)

)(80)

M ha autovalori −a e β(ξ − ξ), la stabilita

dipende dalla mutua posizione di ξ e ξ

29 marzo 2001 79

Per (x, y) = (0, η)

M =

(a− α b

B 0

−β bB b− 2B bB

)=

=

α (aα − bB

)0

−β bB −b

=

=

(α(η − η) 0

−β bB −b

)(81)

M ha autovalori −b e α(η − η), la stabilita

dipende dalla mutua posizione di η e η

29 marzo 2001 80

Per (x, y) = (x, y)

M =

(−Ax −αx−βy −By

)(82)

se indichiamo con D il determinante di M e

con T la traccia di M , avremo che

D = xy(AB − αβ) T = −(Ax+By)(< 0)

(83)

Poiche si ha sempre

t2

4−D = (xA− yB)2 + 4xyαβ > 0

dal momento che x > 0, y > 0, Gli autovalori

di M sono sempre reali e distinti.

Se D > 0 M ha due autovalori con lo stesso

segno di T (soluzione stabile)

Se D < 0 M ha due autovalori reali di segno

opposto. (soluzione instabile)

29 marzo 2001 81

Modello di crescita

di due popolazioni

in cooperazione

29 marzo 2001 82

Cooperazione obbligatoria

1. La popolazione x, isolata, decrescerebbe

secondo la legge

x(t) = −ax(t) (84)

2. La popolazione y, isolata, decrescerebbe

secondo la legge

y(t) = −by(t) (85)

3. Il tasso di crescita della x aumenta pro-

porzionalmente ad y, cosı che

x(t) = (−a+ βy(t))x(t) (86)

4. Il tasso di crescita della y aumenta pro-

porzionalmente ad x, cosı che

y(t) = (−b+ αx(t))y(t) (87)

29 marzo 2001 83

ne consegue che

x(t) = (−a+ βy(t))x(t)

y(t) = (−b+ αx(t))y(t)(88)

Il sistema ammette come soluzioni costanti

(punti critici)

(0,0) , (a

β,b

α)

La stabilita puo essere studiata linearizzando

∇((−a+ βy)x(−b+ αx)y

)=

(−a+ βy βx

αy −b+ βx

)(89)

29 marzo 2001 84

Per x = 0, y = 0, otteniamo

(−a 00 −b

)(90)

La matrice ha autovalori entrambi negativi

La soluzione nulla e asintoticamente stabile.

Per x = aβ, y = b

α si ha

(0 bβαbαβ 0

)(91)

La matrice ha autovalori√ab , −

√ab

La soluzione nulla e instabile.

29 marzo 2001 85

Esempio di cooperazione tra due popolazioni.

29 marzo 2001 86

Cooperazione facoltativa

1. La x, in assenza della y, crescerebbe se-

condo la legge logistica

x(t) = (a− b(x(t))x(t) (92)

2. La y, in assenza della x, crescerebbe se-

condo la legge logistica

y(t) = (c− dy(t))y(t) (93)

3. La presenza della y aumenta il tasso di

crescita della x proporzionalmente ad x,

cosı che

x(t) = (a− b(x(t) + γy(t))x(t) (94)

4. La presenza della x aumenta il tasso di

crescita della y proporzionalmente ad x,

cosı che

y(t) = (c− dy(t) + δx(t))y(t) (95)

29 marzo 2001 87

ne consegue che

x(t) = (a− b(x(t) + γy(t))x(t)

y(t) = (c− dy(t) + δx(t))y(t)(96)

Che ammette come punti critici

E1 = (0,0) , E4 = (A

D,B

D) (97)

E2 = (a

b,0) , E3 = (0,

c

d) (98)

dove

A = cγ + ad (99)

B = cb+ aδ (100)

D = bd− γδ (101)

29 marzo 2001 88

La matrice Jacobiana per E = (x, y) e

M =

(a− 2bx+ γy γx

δy c− 2dy+ δx

)(102)

Il punto E1 = (0,0) e un nodo instabile

M =

(a 00 c

)a, b > 0 (103)

Il punto E2 = (ab ,0) e un punto sella

M =

(−a γab0 c+ δab

)− a < 0 c+ δ

a

b> 0

(104)

29 marzo 2001 89

Il punto E3 = (0, cd) e un punto sella

M =

(a+ γ cd 0γ cd −c

)a+ γ

c

d> 0 − c < 0

(105)

Il punto E4 = (AD,BD) appartiene al primo

quadrante solo se D > 0 e perche il model-

lo abbia senso questa condizione deve essere

soddisfatta. In tal caso il punto e stabile.

M =

(−bAD −γADδBD −dBD

)(106)

e

detM = (bd− γδ)AB = DAB > 0 (107)

mentre la traccia di M e negativa.

29 marzo 2001 90

Esempio di cooperazione tra due popolazioni.

29 marzo 2001 91

Diffusione di un’epidemia

di tipo

SIS o SIR

E un modello simile a quello precedente

descrive la diffusione di una epidemia.

SIS significa

S uscettibili (all’infezione)

I nfetti

S uscettibili (nuovamente all’infezione)

mentre SIR sta ad indicare

S uscettibili

I nfetti

R imossi (colpiti dall’infezione, guariti e im-muni o morti).

29 marzo 2001 92

Sia x il numero di individui infetti

Sia y il numero di individui suscettibili,

L’accrescimento degli infetti e proporzionale

al numero di incontri tra x e y diminuito della

parte di infetti che guariscono.

b e il tasso di guarigione cioe il coefficiente

che indica la frazione guaritiinfetti

a e il tasso di infettivita, cioe la frazioneinfettatiincontri ,

29 marzo 2001 93

avremo

(SIS)

x(t) = ax(t)y(t)− bx(t)

y(t) = −ax(t)y(t) + bx(t)

x(t0) = x0y(t0) = y0

(108)

I guariti (−bx(t)) rientrano nel numero dei

suscettibili (SIS)

mentre nel caso (SIR)

(SIR)

x(t) = ax(t)y(t)− bx(t)

y(t) = −ax(t)y(t)x(t0) = x0y(t0) = y0

(109)

gli individui guariti non sono piu infettabili

29 marzo 2001 94

Nel caso del sistema (SIS),

Il numero degli infetti ed il numero dei su-

scettibili ha somma costante uguale al nu-

mero totale degli individui della popolazione

N

Il sistema puo essere ridotto ad una singola

equazione tenendo conto che

y = N − x

In tal caso si ha

x(t) = ax(t)(N − x(t))− bx(t) (110)

e ci si riduce allo studio della solita equazione

logistica.

29 marzo 2001 95

Nel caso del sistema (SIR) osserviamo che

x(t) = x0e

∫ tt0

(ay(s)−b)ds(111)

y(t) = y0e

∫ tt0−ax(s)ds

(112)

e quindi

x, y > 0

Inoltre sommando membro a membro in (SIR)

si ottiene

x(t) + y(t) = −bx(t) ≤ 0 (113)

t 7→ x(t) + y(t)

e decrescente e quindi limitata da x(t0) +

y(t0) = N percio x, y sono limitate ed esisto-

no per ogni t ≥ t0.

29 marzo 2001 96

Integrando la 113 e ricordando che per la 111

x, y > 0, si ha che

b∫ tt0x(s)ds ≤ N − x(t)− y(t) ≤ N (114)

l’integrale a primo membro esiste in senso

improprio per t→ +∞ ed e convergente.

(si ricordi che x ≥ 0 e che y ≥ 0);

Ne segue che (se x ammette limite all’infini-

to) x(t) → 0 per t→ +∞.

Possiamo rappresentare le soluzioni del siste-

ma (SIR) nel piano (x, y).

29 marzo 2001 97

Se ϕ rappresenta localmente la x in funzio-

ne della y dove x ed y sono le soluzioni del

sistema dato, avremo

x(t) = ϕ(y(t)) = 0

e

x(t) = ϕ′(y(t))y(t)

da cui

ϕ′(y) =axy − bx

−axy= −1 +

b

ay

Se ne ricava

x = −y+b

aln y+ cost (115)

e cio consente di disegnare facilmente le or-

bite del sistema (SIR) nel piano delle fasi.

29 marzo 2001 98

Andamento del numero degli individui infetti x

e del numero degli individui suscettibili y

nel caso di una epidemia di tipo SIR

(piano delle fasi)

29 marzo 2001 99

Andamento del numero degli individui infetti x

e del numero degli individui suscettibili y

nel caso di una epidemia di tipo SIR

29 marzo 2001 100

Modelli

Preda Predatore

di Lotke-Volterra

Descrive l’andamento di due popolazioni una

delle quali si nutre dell’altra.

Si fonda sulle seguenti ipotesi:

x(t) e il numero di individui della popolazione

preda

y(t) e il numero di individui della popolazione

di predatori

29 marzo 2001 101

1. In assenza di predatori x (le pre-

de) seguono un modello malthusiano

(esponenziale)

x(t) = ax(t)

2. Il tasso di crescita a delle prede diminui-

sce in proporzione al numero di predatori

a− αy(t)

3. In assenza di prede i predatori si estin-

guono secondo un modello malthusiano

(esponenziale)

y(t) = −by(t)

4. Il tasso di accrescimento b dei predato-

ri aumenta in maniera proporzionale al

numero di prede −b+ βx(t)

29 marzo 2001 102

la dinamica delle due popolazioni puo esse-

re descritta dal seguente sistema differenziale

(Equazioni di Lotka-Volterra)x(t) = ax(t)− αy(t)x(t)

y(t) = −by(t) + βx(t)y(t)(116)

con i dati inizialix(t0) = x0y(t0) = y0

(117)

ciascuna popolazione cresce proporzionalmen-

te al numero dei suoi individui ed il nume-

ro degli incontri xy con gli individui dell’altra

popolazione aumentano la crescita di una e

diminuiscono la crescita dell’altra

29 marzo 2001 103

1. Il sistema ammette due puntistazionari ( bβ ,

aα), (0,0)

2. Le orbite che giacciono nel primoquadrante sono curve chiuse

3. Le soluzioni del sistema sonoperiodiche

4. Vale per le soluzioni del sistemaun principio di conservazione dellemedie in quanto si vede che, postoSi ha

1

T

∫ T0x(s)ds =

b

β,

1

T

∫ T0y(s)ds =

a

α(118)

29 marzo 2001 104

Il sistema 116 ha soluzioni costanti

(x(t), y(t)) = (ξ, η)

se 0 = aξ − αηξ

0 = −bη+ βηξ(119)

cioe quando

ξ = 0, η = 0

oppure quando

ξ =b

β, η =

a

α

29 marzo 2001 105

Il sistema linearizzato relativo a 116 ha come

matrice dei coefficienti(a− αy −αxβy −b+ βx

)(120)

Per (x, y) = (0,0) la matrice diventa(a 00 −b

)(121)

e si vede che la soluzione nulla e instabile.

Per (x, y) =(bβ ,

aα

)la matrice diventa(0 −ab 0

)(122)

Gli autovalori sono complessi coniugati con

parte reale nulla; il sistema linearizzato e sta-

bile ma non si puo concludere per il sistema

originale.

29 marzo 2001 106

Per ottenere informazioni sulle soluzioni del

sistema possiamo operare in maniera piu di-

retta.

Possiamo cioe tracciare le orbite del sistema,

Ovvero le curve descritte parametricamente

da x = x(t)

y = y(t)(123)

Cio viene fatto considerando un integrale pri-

mo del sistema.

Riscriviamo il sistema nella forma nella formax(t) = ax(t)− αy(t)x(t) = φ(x(t), y(t))

y(t) = −by(t) + βx(t)y(t) = ψ(x(t), y(t))

(124)

29 marzo 2001 107

Cerchiamo soluzioni che possano essere de-

scritte da funzioni y(x) In tal caso si potra

porre x(t) = x

y(t) = y(x(t)) = y(x)

e si potra ricavare che y(t) = y′(x(t))x(t), per

cui sostituendo nella 124

y′(x) =φ(x, y(x))

ψ(x, y(x))=−by(x) + βxy(x)

ax− αxy(x)(125)

29 marzo 2001 108

Separando le variabili si ottiene

y′(x)

(a

y(x)− α

)=(−bx

+ β

)(126)

Integrando

a ln y(x)− αy(x) = −b lnx+ βx+ c (127)

Ora se poniamo

g(y) = a ln y − αy f(x) = −b lnx+ βx

(128)

possiamo riscrivere la 127 come

g(y) = f(x) + c (129)

29 marzo 2001 109

La seguente figura rappresenta l’andamentodella funzione g;

y0 =a

α(130)

g0 = g(y0) = a lna

α− α

a

α= a

(lna

α− 1

)> 0

(131)

Poiche e ragionevole supporre le costanti a, bmolto maggiori, rispettivamente, di α, β pos-siamo assumere che a

α > e

29 marzo 2001 110

La seguente figura rappresenta l’andamentodella funzione f + c;

Si ha x0 = bβ f0 = f(x0) < 0

Se g−1 e l’inversa di una opportuna restrizio-ne di g dalla 129

y(x) = g−1(f(x) + c) (132)

29 marzo 2001 111

Osserviamo (si veda la figura seguente) che129 ha soluzioni se e solo se

[f0 + c, g0] = (−∞, g0] ∩ [f0 + c,+∞) 6= ∅

e che

f(xm) + c = f(xM) + c = g0

g(ym) = g(yM) = f0 + c

29 marzo 2001 112

Inoltre il grafico delle inverse di g (inverse

perche occorre considerare separatamente di-

verse restrizioni) e dato da

29 marzo 2001 113

la funzione f e decrescente per x ∈ [xm, x0]

ed assume tutti i valori compresi tra g0 e f0+

c,

g e crescente ed invertibile sull’intervallo [ym, y0]

ed assume valori compresi tra g0 e f0 + c.

g−1 e definita da [f0+c, g0] a valori in [ym, yM ]

g−1(f + c) risulta definita e decrescente su

[xm, x0]

Analogamente su [x0, xM ] risulta crescente.

29 marzo 2001 114

Tenendo conto anche dell’inversa relativa al-

l’intervallo [y0, yM ] possiamo disegnare le tra-

iettorie del sistema

L’andamento della curva definita dalla 129 si

puo osservare nella seguente figura

29 marzo 2001 115

La funzione U(x, y) = a ln y − αy + b lnx− βx

ottenuta nella 127, si chiama integrale primo

del sistema

Le orbite del sistema costituiscono le sue cur-

ve di livello: infatti si ha, proprio dalla 127

U(x, y(x)) = c (133)

Se

U(x(t), y(t)) = 0 (134)

allora

Ux(x(t), y(t))x(t) + Uy(x(t), y(t))y(t) = 0

(135)

da cui

φ(x(t), y(t))

ψ(x(t), y(t))=x(t)

y(t)= −

Ux(x(t), y(t))

Uy(x(t), y(t))(136)

29 marzo 2001 116

Dalla 136 il gradiente di U e proporzionale al

secondo membro del sistema differenziale in

esame;

Il campo di direzioni definito dal sistema ed

il campo di direzioni definito dal gradiente di

un integrale primo coincidono.

(coincidono le direzioni, non i vettori in quan-

to il fattore di proporzionalita puo non so-

lo non essere 1 ma puo anche non essere

costante.)

Lo studio delle orbite del sistema consente

anche di verificare che

29 marzo 2001 117

1. Le orbite del sistema sono chiuse

2. Le orbite contengono al loro interno il

punto stazionario (x0, y0)

3. Le orbite sono contenute nel primo

quadrante (xm, xM , ym, yM > 0)

4. le orbite non passano per nessun pun-

to stazionario (gli unici punti stazionari

diversi da (x0, y0) sono sugli assi)

Poiche il sistema e autonomo, le traslate di

una soluzione sono ancora soluzioni;

cioe se (x(t), y(t)) risolve il problema allora

anche (x(t+ T ), y(t+ T )) e soluzione.

Per ogni punto del piano delle fasi passa una

ed una sola orbita in virtu di quanto visto

sulle orbite dei sistemi lineari .

29 marzo 2001 118

Le soluzioni del sistema sono periodiche;

infatti poiche la curva e costituita di due

tratti che possono essere rappresentati co-

me funzioni su un intervallo limitato la sua

lunghezza e finita `,

inoltre

x2(t) + y2(t) =

= φ2(x(t), y(t)) + ψ2(x(t), y(t)) >

> min(x− x)2 + (y − y)2 > ε

x > 0 , y > 0

{φ2(x, y)+ψ2(x, y)} > m > 0

(137)

in quanto φ, ψ si annullano contemporanea-

mente solo nei punti (0,0) e (x, y)

29 marzo 2001 119

Pertanto∫+∞0

√x2(t) + y2(t)dt = +∞ ed esi-

ste T > 0 tale che∫ T0

√x2(t) + y2(t)dt = `

Dopo un tempo T la soluzione deve ripassare

per il punto iniziale (x, y)

Per il teorema di unicita, da li’ e costretta a

ripercorrere la stessa orbita,

Percio la soluzione e periodica

29 marzo 2001 120

Si ha inoltrex(t)x(t) = a− αy(t)y(t)y(t) = −b+ βx(t)

(138)

Per la periodicita delle soluzioni si ottiene

0 = lnx(T )

x(0)=∫ T0a− αy(s)ds =

= aT − α∫ T0y(s)ds

(139)

0 = lny(T )

y(0)=∫ T0−b+ βyxs)ds =

= bT − β∫ T0x(s)ds (140)

per cui si puo ricavare la 118

b

β=

1

T

∫ T0x(s)ds

a

α=

1

T

∫ T0y(s)ds

(141)

29 marzo 2001 121

Modelli Preda Predatore

con

Prelievo costante

In un sistema preda predatore consideriamo

l’effetto di un prelievo costante h sulle prede

e k sui predatori.

Siamo ricondotti ax(t) = ax(t)− αy(t)x(t)− hx(t)

y(t) = −by(t) + βx(t)y(t)− ky(t)

x(t0) = x0y(t0) = y0

(142)

Il comportamento delle soluzioni e simile a

quello del sistema senza prelievo se a−h > 0.

Viene modificata la condizione di equilibrio

ed il valore medio

x =b+ k

βy =

a− h

α(143)

29 marzo 2001 122

La legge di Lanchester

Rappresenta l’andamento di due popolazioni

F ed E in lotta

Supponendo nota la capacita di ognuna di

esse di causare danni all’altra.

dF e dE misurano tale capacita

F (t) = −dEE(t)

E(t) = −dFF (t)(144)

Se F = (Friends), E = (Enemies), dE e dFsono i coefficienti di efficienza dei medesimi,

Il modello e dovuto all’ingegnere aeronautico

Lanchester

(Progetto il bombardiere Lanchester e scrisse

diversi lavori su problemi economici.)

29 marzo 2001 123

La soluzione del sistema e facilmente otteni-bile

Tuttavia alcune proprieta possono essere de-dotte dal sistema senza cercarne la soluzioneesplicita.

moltiplicando la prima delle 144 per dFF (t)

e la seconda per dEE(t)

dFF (t)F (t) = −dEdFF (t)E(t)

dEE(t)E(t) = −dFF (t)dEE(t)(145)

sottraendo membro a membro

dFF (t)F (t) = dEE(t)E(t) (146)

da cui

dFF2(t)− dEE

2(t) = cost. (147)

29 marzo 2001 124

Se E0 ed F0 e la consistenza numerica iniziale

dei due eserciti

dFF2 − dEE

2 = dFF20 − dEE

20 (148)

”legge dei quadrati di Lanchester”

Possiamo calcolare la consistenza delle for-

ze dell’esercito vittorioso Ffinale al termine

della guerra supponendo che in tale istante

la consistenza Efinale dell’esercito nemico sia

ridotta a zero.

Avremo

F2finale = F0 −

dEdFE2

0 (149)

Dalla 149 si nota come

l’esercito vittorioso subisce perdite, alla fine

della guerra, proporzionali al quadrato della

consistenza iniziale del nemico.29 marzo 2001 125

La non linearita del risultato puo essere mes-

sa in evidenza dal seguente esempio.

Siano

F = 50000 unita

E1 = 40000 ed E2 = 30000 unita

• L’esercito F affronta in successione prima

E1 e poi E2.

In tal caso, dopo il primo scontro F si

riduce a 30000 unita in quanto

(50000)2 − (40000)2 = (30000)2

mantiene una consistenza pari a quella di

E2 e c’e equilibrio

29 marzo 2001 126

• L’esercito F affronta E = E1 + E2 =

70000.

In tal caso F ha la peggio in quanto

(50000)2 − (70000)2 ≈ −(48990)2

Poiche e facile calcolare esplicitamente la so-

luzione F (t), E(t)) del sistema,

possiamo determinare il tempo necessario al-

l’esercito vincitore F per conquistare la vit-

toria risolvendo

E(t) = 0

29 marzo 2001 127

Il Pendolo

e

Il Pendolo Rovesciato

Il pendolo semplice

Consideriamo una pallina P di massa m at-

taccata all’estremita di un’asta di massa tra-

scurabile imperniata nell’altra estremita.

29 marzo 2001 128

Se

1. ` e la lunghezza dell’asta

2. θ(t) e l’angolo che l’asta forma con la

verticale nell’istante t

3. s(t) = `θ(t) e la distanza, sull’arco di cir-

conferenza descritto dalla pallina, dalla

posizione allineata con la verticale

29 marzo 2001 129

avremo che

s(t) = `θ(t) (150)

s(t) = `θ(t) (151)

s(t) = `θ(t) (152)

La pallina si suppone sottoposta

• alla forza peso

• ad una forza che si oppone alla direzio-

ne del moto proporzionale alla velocita

(resistenza del mezzo).

Pertanto l’equazione che governa il sistema

sara data da

m`θ(t) = −k`θ(t)−mg sin θ(t) (153)

29 marzo 2001 130

A meno di ridenominare le costanti in gioco

si ottiene

θ(t) = −Rθ(t)− P sin θ(t) (154)

Posto x(t) = θ(t)

y(t) = θ(t)(155)

otteniamo che il sistema e governato dalle

equazionix(t) = y(t)

y(t) = −Ry(t)− P sinx(t)(156)

29 marzo 2001 131

Ci sono due punti di equilibrio

(x, y) = (0,0) e (x, y) = (π,0)

La loro stabilita si puo studiare usando il

sistema linearizzato che:

• Per (x, y) = (0,0) ex(t) = y(t)

y(t) = −Ry(t)− Px(t)(157)

• Per (x, y) = (0, π) e(x(t)− π)′ = y(t)

y(t) = −Ry(t) + P (x(t)− π)(158)

Lo studio degli autovalori delle matrici dei

coefficienti mostra che per il primo sistema

(0,0) e un punto di equilibrio stabile mentre

per il secondo sistema (π,0) e un punto di

equilibrio instabile.

29 marzo 2001 132

Il pendolo rovesciato

Poniamoci ora nella situazione dell’esempio

precedente ma consideriamo la variabile θ mi-

surata a partire dalla posizione verticale.

Il pendolo rovesciato

In questo nuovo riferimento avremo che la

posizione di equilibrio (x, y) = (0,0) e in-

stabile (corrisponde alla precedente (π,0)) e

possiamo porci il seguente problema

29 marzo 2001 133

E possibile, applicando una op-

portuna forza sul perno dell’asta,

mantenere la pallina in equilibrio?

Applichiamo sul perno una forza v = −mw(t)

rivolta in modo da contrastare la caduta della

pallina; la forza applicata si puo rappresenta-

re agente in P .

L’equazione che governa il sistema e data da

m`θ(t) = −k`θ(t) +mg sin θ(t)−mw(t) cos θ(t)

(159)

da cui

θ(t) = −k

m`θ(t) +

g

`sin θ(t)−

w(t)

`cos θ(t)

(160)

29 marzo 2001 134

Posto x(t) = θ(t)

y(t) = θ(t)(161)

Il sistema e governato dalle equazionix(t) = y(t)

y(t) = − km`y(t) + g

` sinx(t)−w(t)` cosx(t)

(162)

e ancora, a meno di ridenominare le costanti,x(t) = y(t)

y(t) = −Ry(t) + P sinx(t)− w(t)` cosx(t)

(163)

29 marzo 2001 135

Se supponiamo che la forza riequilibrante il

sistema sia lineare in θ e θ,

w(t)

`= aθ(t) + bθ(t) = ax(t) + by(t)

Il sistema diventax(t) = y(t)y(t) =−Ry(t) + P sinx(t)−

− (ax(t) + by(t)) cosx(t)(164)

e linearizzando in (0,0y)x(t) = y(t)

y(t) = (−b−R)y(t) + (−a+ P )x(t)

(165)

29 marzo 2001 136

Possiamo scegliere a e b in modo che il siste-

ma sia stabile

Posto b+R = µ

P − a = η(166)

Il sistema si scrivex(t) = y(t)

y(t) = ηx(t)− µy(t)(167)

Gli autovalori del sistema sono dati da

−µ±√µ2 + 4η

2ed e facile determinare µ e η in modo che il

sistema sia stabile.

Ad esempio µ > 0 e η < 0, ma non solo.

29 marzo 2001 137

Modelli differenziali

collegati alla

conservazione della massa

In un sistema isolato il bilancio tra ”la quan-

tita di materia che entra, quella che esce e

quella che che rimane” deve essere in pareg-

gio.

Si possono costruire

1. Modelli differenziali ordinari

2. Modelli differenziali alle derivate parziali

29 marzo 2001 138

Modelli differenziali ordinari

Capacita di un lago

L’acqua entra in un lago con un flusso

costante di K litri al minuto e dal lago

ne evapora una quantita proporziona-

le secondo una costante H a v2/3 do-

ve con v si indica il volume di acqua

presente nel lago

29 marzo 2001 139

Sia v(t) il volume di acqua nel lago all’istante

t,

lim∆t→0

v(t+ ∆t)− v(t)

∆t= v(t) (168)

In assenza di evaporazione

v(t) = K (169)

onde

v(t) = v(0) +Kt

Se c’e evaporazione (perdita di volume pro-

porzionale a v2/3(t)) pari a

−Hv2/3(t)

avremo

v(t) = K −Hv2/3(t) (170)

29 marzo 2001 140

l’equazione 170 ammette una soluzione co-stante

v(t) =(K

H

)3/2

E soluzione di equilibrio?

Andamento della quantita di acqua nel lago

29 marzo 2001 141

Poiche se v(t) <(KH

)3/2allora v(t) > 0

quindi v(t) e crescente,

Mentre se v(t) >(KH

)3/2allora v(t) < 0

e quindi v(t) e decrescente;

La soluzione costante puo essere considerata

la soluzione di regime.

La soluzione costante e asintoticamente sta-

bile.

29 marzo 2001 142

Grado di inquinamento di un lago

Un secondo esempio si puo costruire consi-

derando la seguente situazione

Un inquinante entra in un lago con flusso

costante σ, p(t)e la massa di inquinante al

tempo t; esso viene metabolizzato dai batte-

ri presenti nel lago in quantita proporziona-

le alla sua massa secondo una costante k In

questo processo viene consumata una quan-

tita dell’ossigeno disciolto nelle acque del la-

go pari alla massa di inquinante decomposto.

Il livello o(t) di ossigeno nel lago tuttavia e

reintegrato attraverso il contatto tra la su-

perficie dell’acqua e l’aria in maniera propor-

zionale alla differenza om − o(t) tra il valore

di saturazione dell’ossigeno om e la quantita

di ossigeno presente o(t)

29 marzo 2001 143

Le equazioni che descrivono questo fenome-

no possono pertanto essere scritte nella se-

guente manierap(t) = σ − kp(t)

o(t) = −kp(t) + h(om − o(t))(171)

Indichiamo le soluzioni di equilibrio della 171

come

p =σ

ko = om −

σ

h(172)

Possiamo linearizzare il sistema in (p, o) con

semplici operazioni algebriche, (e gia lineare

ma non omogeneo) ottenendop = −k(p− σ

k) = −k(p− p)

o =− k(p−σ

k+σ

k) +

+ h(om − (o− (om −σ

h) + (om −

σ

h)))

(173)

29 marzo 2001 144

da cuip = −k(p− p)o =− k(p− p)− σ+

+h(om − (o− o)− om +σ

h)

(174)

ed ancorap = −k(p− p)

o = −k(p− p)− h(o− o)(175)

la matrice dei coefficienti del sistema e

A =

(−k 0−k −h

)(176)

ha autovalori −k, −h, negativi e garantiscono

la stabilita.

29 marzo 2001 145

Andamento dell’inquinante e dell’ossigeno

(piano delle fasi)

29 marzo 2001 146

Andamento dell’inquinante e dell’ossigeno

29 marzo 2001 147

Modelli differenziali

alle derivate parziale

Vogliamo studiare i dati relativi alla distribu-

zione di una massa che si muove lungo l’asse

x nel tempo t.

• ρ(x, t) e la densita della sostanza che in-

tendiamo studiare,

• x e la coordinata relativa all’asse su cui

avviene il movimento

• t rappresenta la variabile di tempo.

29 marzo 2001 148

Le dimensioni di ρ sono unita di massaunita di lunghezza

la massa compresa tra x e x+ ∆x e

m(t, x, x+ ∆x) =∫ x+∆x

xρ(t, s)ds

φ(x, t) e il flusso attraverso x al tempo t,

cioe

φ(x, t) e la quantita di massa che transita per

il punto x nell’istante t;

le dimensioni di φ sono unita di massaunita di tempo .

29 marzo 2001 149

Consideriamo un tratto tra x e x+∆x e due

istanti successivi t e t+ ∆t.

Sia m(t, x) la massa presente all’istante t nel

punto x

La massa in [x, x+ ∆x] all’istante t e∫ x+∆x

xρ(t, s)ds = ρ(t, x)∆x x ≤ x ≤ x+ ∆x

(177)

La massa che transita attraverso il punto x

nel tempo [t, t+ ∆t] e∫ t+∆t

tφ(s, x)ds = φ(t, x)∆t t ≤ t ≤ t+ ∆t

(178)

29 marzo 2001 150

Consideriamo un intervallo di lunghezza I =[x, x + ∆x] ed un intervallo di tempo T =[t, t+ ∆t]

1. al tempo t

2. e al tempo t+ ∆t

Avremo che

1. la massa in I all’istante t e ρ(t, x)∆x

2. la massa in I all’istante t + ∆t e ρ(t +∆t, x)∆x

3. la massa che transita in x nel tempo trat e t+ ∆t e φ(t, x)∆t

4. la massa che transita in x+∆x nel tempotra t e t+ ∆t e φ(t, x)∆t

29 marzo 2001 151

Uguagliando la variazione della massa conte-

nuta in I nel tempo T con la differenza di

quanto e entrato in I, durante il tempo T(ρ(t+ ∆t, x)− ρ(t, x)

)∆x =

=(φ(t, x)− φ(t, x+ ∆x)

)∆t (179)

Dividendo per ∆t∆x

ρ(t+ ∆t, x)− ρ(t, x)

∆t=φ(t, x)− φ(t, x+ ∆x)

∆x(180)

e passando al limite per ∆x → 0 e ∆t → 0

poiche

x ≤ x, x ≤ x+ ∆x t ≤ t, t ≤ t+ ∆t

(181)

si ha

∂ρ(x, t)

∂t= −

∂φ(x, t)

∂x(182)

29 marzo 2001 152

Alternativamente poiche il flusso e la derivatadella massa rispetto al tempo

φ(t, x) =d

dtm(t, x) (183)

e

φ(t, x+ ∆x) =d

dtm(t, x+ ∆x) (184)

quindi

φ(t, x)− φ(t, x+ ∆x)

∆x=

=1

∆x

d

dt(m(t, x)−m(t, x+ ∆x)) =

= −1

∆x

d

dt

∫ x+∆x

xρ(t, s)ds (185)

e scambiando derivata ed integrale

φ(t, x)− φ(t, x+ ∆x)

∆x=

−1

∆x

∫ x+∆x

x

(∂

∂tρ(t, s)

)ds (186)

Passando al limite per ∆x → 0 si ottienel’equazione 182. ∆x→ 0.

29 marzo 2001 153

Si puo tenere conto di apporto o sottrazione

di massa

mediante un termine k(x, t)

∂ρ(x, t)

∂t= −

∂φ(x, t)

∂x+ k(x, t) (187)

L’equazione 187 contiene troppe incognite

(ρ e φ) e quindi occorre introdurre qualche

ulteriore condizione;

29 marzo 2001 154

Osserviamo che la 187 e stata ricavata dalla

180 passando al limite per ∆x→ 0 e ∆t→ 0.

Tuttavia si e trascurato il fatto che ∆t e ∆x

possono non essere indipendenti.

Tale omissione si corregge con ulteriori con-

dizioni

Ad esempio possiamo supporre che

φ = φ(ρ) (188)

cosı che

∂φ

∂x=

(dφ

dρ

)(∂ρ

∂x

)(189)

l’equazione 187 diventa allora

∂ρ(x, t)

∂t= −

(∂φ(ρ(x, t))

∂ρ

)(∂ρ(x, t)

∂x

)+ k(x, t)

(190)

29 marzo 2001 155

Possiamo poi considerare due casi partico-

larmente significativi di dipendenza di φ da

ρ

se φ(ρ) = v(x, t)ρ il fenomeno e di advezione.

v(x, t) ha dimensione unita di lunghezzaunita di tempo (di-

mensioni di una velocita)

e la velocita con cui la materia si muove lungo

l’asse x;

infatti

v = φρ = unita di massa

unita di tempounita di lunghezza

unita di massaAd esempio v(x, t) = c

quando

φ(ρ) = −ν∂ρ

∂x(x, t)

si parla di diffusione

E caratteristica della propagazione del calore.

Il calore tende a fluire, proporzionalmente al

gradiente di temperatura (la densita di calo-

re), dalla temperatura piu alta verso la piu

bassa.

29 marzo 2001 156

La tipica equazione di advezione prende la

forma

∂ρ(x, t))

∂t= −

∂(ρ(x, t)v(x, t))

∂x+ k(x, t) (191)

La tipica equazione di diffusione e

∂ρ(x, t)

∂t= ν

∂2ρ(x, t)

∂x2+ k(x, t) (192)

Piu in generale gli effetti di advezione e di

diffusione possono sovrapporsi ed in tal caso

l’equazione diventa

∂ρ(x, t)

∂t= ν

∂2ρ(x, t)

∂x2−∂(ρ(x, t)v(x, t)

∂x+ k(x, t)

(193)

29 marzo 2001 157

Modello di trasporto

di

un inquinante in un fiume

Consideriamo la seguente situazione

Un inquinante organico e mescolato alle ac-

que di un fiume, che scorre con velocita co-

stante c lungo l’asse x; vogliamo conosce-

re la concentrazione ρ(x, t), che supponia-

mo omogenea in ciascuna sezione del fiume,

dell’inquinante nota la distribuzione iniziale

ρ(x,0) dell’inquinante stesso, tenendo conto

che l’inquinante viene degradato dall’azione

batterica proporzionalmente alla concentra-

zione di inquinante e trascurando i fenomeni

diffusivi.

29 marzo 2001 158

La situazione appena descritta e un esempio

di advezione: la velocita con cui si muove la

massa e c.

Puo essere modellizzata utilizzando l’equa-

zione

∂

∂tρ(x, t) = −c

∂

∂xρ(x, t)− µρ(x, t) (194)

µρ(x, t) descrive la metabolizzazione dell’in-

quinante

29 marzo 2001 159

Riscriviamo la 194 nella forma

∂

∂tρ(x, t) + c

∂

∂xρ(x, t) = −µρ(x, t) (195)

Il primo membro e la derivata rispetto a t

della funzione

R(x, t) = ρ(x+ ct, t) (196)

e quindi

d

dtR(x, t) =

d

dtρ(x+ ct, t) =

=∂

∂tρ(x+ ct, t) + c

∂

∂xρ(x+ ct, t) =

= −µρ(x+ ct, t) = −µR(x, t) (197)

29 marzo 2001 160

Integrando rispetto a t l’equazione differen-

ziale lineare che ne risulta

ρ(x+ ct, t) = h(x)e−µt (198)

Imponendo le condizioni iniziali ρ(x,0) = ρ0(x)

(e nota la densita iniziale in ogni luogo del

fiume)

ρ0(x) = ρ(x,0) = h(x) (199)

da cui

ρ(x+ ct, t) = ρ0(x)e−µt (200)

e

ρ(x, t) = ρ0(x− ct)e−µt (201)

29 marzo 2001 161

Consideriamo la situazione in cui

Si verifichi in una locazione x = 0 che possia-

mo supporre coincidente con l’origine una im-

missione di inquinante con flusso costante γ.

(Supponendo che il fiume sia in precedenza

pulito)

Allo scopo, posto

γ(t) =

0 t < 0

γ0 t ≥ 0(202)

possiamo imporre nella 200 che

γ(t) = ρ(0, t) = ρ0(−ct)e−µt (203)

29 marzo 2001 162

Da cui, posto x = −ct si ricava t = −xc

e quindi

γ(−x

c) = ρ0(x)e

µx/c ρ0(x) = γ(−x

c)e−µx/c

(204)

Se ne conclude

ρ(x, t) = ρ0(x− ct)e−µt = γ(ct− x

c)e−µ(x−ct)/ce−µt

(205)

Ed infine

ρ(x, t) = γ(t−x

c)e−µx/c (206)

In questa espressione possiamo osservare che

si evidenzia il fatto che per t− xc < 0 cioe per

x > ct l’inquinante ha densita ρ nulla. (L’in-

quinante, in quel momento non ha ancora

raggiunto quel punto del fiume).

29 marzo 2001 163

Per studiare il livello δ dell’ossigeno disciolto

nel fiume, possiamo tenere conto che anche

δ deve soddisfare una equazione di trasporto

del tipo

∂

∂tδ(x, t) + c

∂

∂xδ(x, t) =

= −µρ(x, t) + µ1(δm − δ(x, t)) (207)

dove −µρ(x, t) rende conto del fatto che il

livello di ossigeno diminuisce in maniera pro-

porzionale all’inquinante metabolizzato,

e µ1δ(x, t) e il reintegro mediante scambio

con l’atmosfera, proporzionale alla differen-

za tra il livello δm di saturazione ed il livello

attuale δ(x, t).

29 marzo 2001 164

L’equazione 207 puo essere integrata come

le precedenti.

Avremo allora che, posto D(t) = δ(x+ ct, t)

dalla 207 si ricava

D′(t) =∂

∂tδ(x+ ct, t) + c

∂

∂xδ(x+ ct, t) =

= −µρ(x+ ct, t) + µ1(δm − δ(x+ ct, t)) =

= −µ1D(t) + µ1δm − µρ0(x)e−µt (208)

ed infine

D′(t) = −µ1D(t) + µ1δm − µρ0(x)e−µt

(209)

Si tratta di una equazione lineare la cui solu-

zione e

D(t) = k(x)e−µ1t + δm +µρ0(x)

µ− µ1e−µt (210)

29 marzo 2001 165

Tenuto conto che il fiume all’istante iniziale

e pulito, possiamo supporre che δ(0, t) = δm,

quindi si puo ricavare k(x) e la 210 diventa

D(t) =µρ0(x)

µ− µ1(e−µt − e−µ1t) + δm (211)

Poiche D(t) = δ(x+ ct, t) possiamo ricavare

che

δ(x+ ct, t) =µρ0(x)

µ− µ1(e−µt − e−µ1t) + δm

(212)

da cui

δ(x+ ct, t) =

=µγ(−xc)e

−µx/c

µ− µ1(e−µt − e−µ1t) + δm (213)

E

δ(x, t) =µγ(t− x

c)e−µ(x−ct)/c

µ− µ1(e−µt − e−µ1t) + δm

(214)

29 marzo 2001 166

Infine

δ(x, t) =µγ(t− x

c)e−µx/c

µ− µ1(1− e(µ−µ1)t) + δm

(215)

Andamento della concentrazione di inquinante

29 marzo 2001 167

Modelli differenziali

di

Traffico autostradale

Il traffico lungo un’autostrada puo essere de-

scritto mediante un modello di trasporto in

cui il flusso dipende dalla densita delle auto.

La natura del modello e discreta

Se consideriamo la situazione ”da molto di-

stante” possiamo ritenere il flusso delle au-

to come il flusso di una massa continua con

densita ρ.

Variabili in gioco

x, descrive l’asse su cui si trova l’autostrada,

t. rappresenta il tempo.

29 marzo 2001 168

L’equazione che governa il flusso delle auto

e di tipo advettivo,di trasporto,

Non esistono ne’ ”sorgenti” ne’ ”pozzi”

Occorre ipotizzare il comportamento del flus-

so φ rispetto alla densita ρ delle auto.

Possiamo fare riferimento a dati che proven-

gono dall’osservazione del traffico in base ai

quali.

1. La velocita delle auto e compresa tra 0

e vm, quest’ultimo valore essendo con-

seguenza dei limiti di velocita o delle

condizioni stradali.

2. La velocita delle auto dipende dalla

densita e se ρm e la massima densita

possibile, supponiamo che

v(ρ) = vm

(1−

ρ

ρm

)(216)

29 marzo 2001 169

Dalla 216 si ricava

• se ρ = 0 allora v = vm

• se ρ = ρm allora v = 0.

L’equazione che ne risulta sara del tipo 190

o piu precisamente del tipo 191

Tenuto conto che

φ(ρ) = v(ρ)ρ = vm

(ρ−

ρ2

ρm

)(217)

da cui

φ′(ρ) =d

dρφ(ρ) = vm

(1−

2ρ

ρm

)(218)

29 marzo 2001 170

L’equazione che regola lo scorrere del traffico

autostradale potra essere scritta nella forma

∂ρ(x, t)

∂t= −vm

(1−

2ρ(x, t)

ρm

)∂ρ(x, t)

∂x(219)

o piu brevemente

ρt = −vm(1−

2ρ

ρm

)ρx = −φ′(ρ)ρx (220)

Cercheremo di studiare la soluzione dell’e-

quazione individuandone le curve di livello;

cercheremo cioe di trovare le curve descritte

dalle equazioni (x(t), t) sulle quali risulta

ρ(x(t), t) = costante = ρ0 (221)

29 marzo 2001 171

Se

ρ(x(t), t) = costante = ρ0 (222)

avremo che

φ′(ρ) = φ′(ρ0) = φ′0

e esso pure costante

Pertanto possiamo riscrivere la 220 come se-

gue

ρt = −φ′0ρx (223)

D’altro canto derivando la 221, si ha

ρt + x(t)ρx = 0 (224)

29 marzo 2001 172

Da 224 e 223 possiamo dedurre che

x(t) = φ′o (225)

da cui, ricordiamo che φ′0 = φ′(ρ0)

x(t) = φ′0t+ x0 (226)

Le curve di livello (x(t), t) della soluzione ρ(x, t)

sono rette di equazione x = φ′0t+ x0.

Se ρ(x,0) = ρ0(x)

si ha

ρ(φ′0t+ x0, t) = ρ(x0,0) = ρ0(x0) (227)

E se supponiamo nota la densita iniziale ρ0(x)

ρ(x, t) = ρ(x− φ′(ρ0)t) = ρ0(x− φ′0t) (228)

29 marzo 2001 173

Esaminiamo qualche esempio

Le auto sono incolonnate ad un semaforo pri-

ma del quale la densita e massima e dopo il

quale la densita e nulla.

Indichiamo con 0 l’istante in cui il semaforo

diventa verde

Il flusso di traffico successivo puo essere de-

scritto dall’equazione 219 con la condizione

iniziale

ρ(x,0) = ρ0(x) =

ρm x < 0

0 x > 0(229)

29 marzo 2001 174

In corrispondenza di tali dati iniziali avremo

che

φ′0 = φ′0(ρ(x,0)) = φ′0(ρ0(x)) =

=

vm x > 0 (ρ0 = 0)

−vm x < 0 (ρ0 = ρm)(230)

Le rette su cui risulta costante la densita,

caratteristiche dell’equazione, sono

x =

vmt+ x0 x0 > 0

−vmt+ x0 x0 < 0(231)

Su ognuna di tali rette la densita e ρ0 ed il

flusso corrispondente e

φ′0 = φ′(ρ(x,0)) = φ′(ρ0)

29 marzo 2001 175

Tali rette coprono solo una parte del semi-

piano t > 0

Non forniscono nessuna informazione su quan-

to accade nella zona di tale semipiano che e

compresa tra le rette x = ±vmt

29 marzo 2001 176

Grafico della funzione densita ρ(x, t)

29 marzo 2001 177

La zona non coperta corrisponde al caso x0 =

0

Non e possibile dare un valore della densita

per x0 = 0

A destra di zero la densita e nulla

A sinistra e massima;

La densita presenta in zero una discontinuita

di tipo ”salto”.

In tale zona ci interessa definire una soluzione

ρ(x, t) in grado di raccordare il valore ρm con

il valore 0.

29 marzo 2001 178

Per x0 = 0 le curve di livello della densita ρ

sono

x = φ′(ρ(0,0))t = φ′(ρ0(0))t (232)

ma la densita iniziale ρ0(0) non e definita,

possiamo soltanto affermare che

ρ0(0) ∈ [0, ρm]

Quindi

φ′(ρ0(0)) = vm

(1−

2ρ0(0)

ρm

)∈ [−vm, vm]

Pertanto per x0 = 0 possiamo considera-

re non una ma infinite rette sulle quali ρ e

costante.

Se in t = 0 assumiamo un valore della densita

ρ, tale valore si manterra costante sulla retta

x = φ′(ρ)t = vm

(1−

2ρ

ρm

)t (233)

29 marzo 2001 179

Ricavando ρ si ottiene

ρ =1

2ρm

(1−

x

vmt

)(234)

Osserviamo che per t = 0 la 234 non e defi-

nita.

Grafico e livelli della 234

29 marzo 2001 180

La 234 soddisfa l’equazione 219; infatti si ha

ρt =x

2vmt2, ρx =

−1

2vmt(235)

Per cui, trascurando la costante 12ρm

x

vmt2+ vm

(1−

2ρ

ρm

)(−1

vmt

)=

=x

vmt2−

1

t+

2ρ

ρmt=

=x

vmt2−

1

t+ρm

(1− x

vmt

)ρmt

= 0 (236)

La soluzione cosı definita ha linee di livelloche sono rette per l’origine

Esse corrispondono ai diversi valori di densitaassunti nell’origine;

Tali valori vanno da

• ρm, in corrispondenza del quale x = −vmt,

• a 0 caso in cui x = vmt

29 marzo 2001 181

Infatti

ρ =1

2ρm

(1−

x

vmt

)= k ⇐⇒

x = vm

(1−

2k

ρm

)(237)

Le rette (al variare di k ∈ [0, ρm] hanno pen-

denza crescente da −vm a vm,

Coprono la zona scoperta dai precedenti cal-

coli

Sono le caratteristiche che passano per l’ori-

gine.

I punti della retta x = −vmt caratterizzano i

tempi ed i luoghi in cui inizia il movimento di

un’auto in coda al semaforo.

Al tempo t iniziano a muoversi le auto che

per t = 0 si trovano alla posizione x = −vmtmentre al tempo t l’auto che si trova in x = 0

per t = 0 avra raggiunto la posizione x = vmt.

29 marzo 2001 182

La velocita con cui le auto si muovono nella

zona in esame, sara

v = vm

(1−

ρ

ρm

)(238)

e quindi dalla soluzione 234 trovata per ρ,

v = vm

(1−

1

2

(1−

x

vmt

))(239)

per cui

v = vm

(1

2+

x

2vmt

)(240)

e

v =vm

2+

x

2t(241)

29 marzo 2001 183

Se x(t) e la posizione di un’auto avremo che

x(t) = v.

Ne viene che

x(t) =vm

2+

x

2t(242)

Si e in questo modo trovata una equazio-

ne differenziale che definisce il movimento

dell’auto alla partenza dopo il verde.

Una condizione iniziale puo essere dedotta te-

nendo conto che l’auto comincia a muoversi

da x0 al tempo t0 = −x0/vm;

Ne viene x(t) = x(t)2t + vm

2x(−x0

vm) = x0

(243)

29 marzo 2001 184

L’integrale generale dell’equazione e

x(t) = C1√t+ vmt (244)

e imponendo che

x0 = x(−x0vm

) = C1

√−x0vm

+ vm

(−x0vm

)(245)

si ricava

C1 = −2√−x0vm (246)

per cui la cui soluzione del problema di Cau-

chy e

x(t) = −2√−x0vm

√t+ vmt =√

t√vm

(√vm√t− 2

√−x0

)(247)

29 marzo 2001 185

Si calcola che l’auto raggiungera il semaforo

al tempo t tale che x(t) = 0; si ricava

t = −4x0vm

La figura seguente mostra come si muove

un’auto in coda al semaforo dopo che il se-

maforo e diventato verde.

29 marzo 2001 186

Una situazione opposta si verifica quando la

densita di traffico aumenta

In tal caso il metodo delle caratteristiche pre-

senta inconvenienti

Si rendono necessari degli aggiustamenti.

La densita iniziale sia data da

ρ(x,0) = ρ0(x) =

ρm4 x < 0

ρm x > 0(248)

Avremo che

φ′0 = φ′0(ρ(x,0)) = φ′0(ρ0(x)) =

=

vm2 x < 0 (ρ0 = ρm4 )

−vm x > 0 (ρ0 = ρm)(249)

29 marzo 2001 187

Le caratteristiche su cui risulta costante la

densita dell’equazione, saranno date da

x =

vm2 t+ x0 x0 < 0

−vmt+ x0 x0 > 0(250)

come si vede nella figura seguente.

29 marzo 2001 188

Grafico della funzione densita ρ(x, t)

29 marzo 2001 189

Dalla figura si vede che la situazione non e

affatto chiara

Esiste una zona del piano in cui le caratteri-

stiche si sovrappongono

Cio causa la mancanza di unicita della solu-

zione, o meglio la sua indeterminatezza.

Questo e dovuto alla discontinuita del dato

iniziale

Dobbiamo operare una scelta tra le due so-

luzioni.

Poiche la densita e discontinua la 182 puo

causare problemi.

29 marzo 2001 190

Ad esempio puo non essere lecito lo scambio

tra derivata ed integrale in 185.

Pertanto ripartiamo dalla 185 stessa.

φ(t, x)− φ(t, x+ ∆x) = −d

dt

∫ x+∆x

xρ(t, s)ds

(251)

Se σ(t) ∈ [x, x + ∆x] e un punto in cui si

verifica la discontinuita la 185 applicata agli

intervalli [x, σ(t)] e [σ(t), x+ ∆x] assicura

d

dt

∫ σ(t)x

ρ(s, t)ds+d

dt

∫ x+∆x

σ(t)ρ(s, t)ds =

= φ(σ+(t), t)−φ(x+∆x, t)−φ(σ−(t), t)+φ(x, t)

(252)

29 marzo 2001 191

D’altro canto per le solite regole di deriva-

zione si ha

d

dt

∫ σ(t)x

ρ(s, t)ds+d

dt

∫ x+∆x

σ(t)ρ(s, t)ds =

σ(t)ρ(σ−(t), t)+∫ σ(t)x

∂ρ(s, t)

∂tds−σ(t)ρ(σ+(t), t)+

+∫ x+∆x

σ(t)

∂ρ(s, t)

∂tds (253)

29 marzo 2001 192

Quindi se x→ σ(t)− e x+ ∆x→ σ(t)+ otte-

niamo

σ(t)ρ(σ−(t), t)− σ(t)ρ(σ+(t), t) =

φ(σ+(t), t)− φ(σ−(t), t) (254)

o piu brevemente

σ(t)ρ− − σ(t)ρ+ = φ+ − φ− (255)

ovvero

σ(t) =φ+ − φ−ρ+ − ρ−

=ρ+v+ − ρ−v−ρ+ − ρ−

(256)

Nel nostro caso si verifica subito che

ρ+ = ρm , ρ− =ρm

4

v+ = 0 , v− =3vm4

29 marzo 2001 193

Per cui

σ(t) = −vm

4e σ(t) = −

vm

4t (257)

La 257 rappresenta l’equazione della curva

lungo la quale si passa da velocita v− = 3vm4

a velocita v+ = 0 ed ivi si riscontra una bru-

sca interruzione del traffico; quindi σ(t) indi-

vidua il punto in cui in traffico subisce uno

shock e ne descrive l’andamento nel tempo;

solitamente σ(t) viene indicata come ”shock

wave” (onda d’urto).

29 marzo 2001 194

Curve

di

Inseguimento

Sia T un bersaglio che si muove su una curva

nota Γ a velocita costante.

Cerchiamo di determinare una curva che deve

percorrere un punto P che insegue T a velo-

cita costante sempre diretta verso il punto

T .

La curva che descritta da P si chiama curva

di inseguimento

Consideriamo alcuni esempi in cui il problema

puo essere risolto.

Il bersaglio T si muove sulla retta x = a a

velocita costante w

29 marzo 2001 195

Γ puo essere parametrizzata mediante le

γ(t) =

x(t) = a

y(t) = wtt ≥ 0 (258)

Siano (x(t), y(t)) le equazioni parametriche

della curva di inseguimento avremo che (x(t), y(t))

deve essere parallelo alla direzione (a−x(t), wt−y(t)) dovra aversi

x(t) = K(a− x(t))

y(t) = K(wt− y(t))

(x(t))2 + (y(t))2 = v2(259)

29 marzo 2001 196

Se supponiamo che la soluzione sia una fun-

zione y(x) potremo affermare che y(t) = y(x(t))

e

y(t) = y′(x(t))x(t) (260)

da cui, tenendo conto della 259

wt− y(t) = y′(x(t))(a− x(t)) (261)

e derivando rispetto t

w − y(t) = y′′(x(t))x(t)(a− x(t))− y′(x(t))x(t)(262)

e, per la 260

w = y′′(x(t))x(t)(a− x(t)) (263)

29 marzo 2001 197

D’altro canto dalla seconda delle 259 possia-

mo ricavare che(dx

dt

)2+(dy

dt

)2=(dx

dt

)2[1 +

(dy

dx

)2]= v2

(264)

Ricavando x(t) dalla 263 e sostituendo nella

264 otteniamo

1 +(y′(x)

)2=(v

w

)2(a− x)2

(y′′)2

(265)

Ora, se poniamo p = y′ e chiamiamo k =vw possiamo ricondurci ad una equazione a

variabili separabili della forma

1 + p2 = k2(a− x)2(p′)2

(266)

29 marzo 2001 198

La curva di inseguimento nel caso in cui

il bersaglio si muova in linea retta

29 marzo 2001 199

Separando le variabili ed integrando ottenia-

mo

p′√1 + p2

=1

k(a− x)(267)

e

(sinh)−1(p) = −1

kln(a− x) + c (268)

da cui

p = sinh(ln(γ(a− x))−1k ) (269)

ed infine

p =1

2

[(γ(a− x))

−1k − (γ(a− x))

1k

](270)

29 marzo 2001 200

La costante γ puo essere ricavata imponen-

do che p(0) = y′(0)) = 0 (Inseguitore inizial-

mente fermo)

Si ottiene che γ = 1a

e possiamo scrivere che

p =1

2

[((1−

x

a))

−1k − ((1−

x

a))

1k

](271)

Integrando nuovamente si puo ottenere la y,

Se k 6= 1

y =1

2

[ka

1− k((1−

x

a))1−

1k −

ka

1 + k((1−

x

a))1+1

k

]+d

(272)

29 marzo 2001 201

Imponendo y(0) = 0 possiamo ricavare d =

− ka1−k2 e ricavare

y =ka

k2 − 1

1 +1

2

(k − 1)(1−

x

a

)1+1k −

−(k+ 1)(1−

x

a

)1−1k

(273)

Se k = 1

p =1

2

[1

γ(a− x)− γ(a− x)

](274)

Imponendo p(0) = y′(0)) = 0 si ottiene γ = 1a

e

p =1

2

[1

1− xa

− (1−x

a)

](275)

Integrando

y =a

4

((1−

x

a

)2− ln

(1−

x

a

)2− 1

)(276)

29 marzo 2001 202

Se k > 1 ( v > w) la 273 permette di calco-

lare che per x→ a y(x) → kak2−1

che pertanto

rappresenta il punto in cui il bersaglio viene

raggiunto dall’inseguitore.

Curva di pedinamento

Se modifichiamo il problema richiedendo che

si mantenga costante la distanza tra insegui-

tore P e bersaglio T ,

(wt− y(t))2 + (a− x(t))2 = a2 (277)

otteniamo sostituendo la 261 nella 277(y′(x)

)2=

a2

(a− x)2− 1 (278)

29 marzo 2001 203

la soluzione e data da

y = −√a2 − (a− x)2 +

+ a ln

a+

√a+

√a2 − (a− x)2

a− x

+ c (279)

Per la 278 la distanza tra inseguitore ed in-

seguito rimane costante l’inseguitore segue,

ma non vuole intercettare, cioe pedina, il

bersaglio.

La stessa equazione vale se l’inseguito traina

l’inseguitore;

la curva descritta dalla 279 si chiama trattri-

ce.

29 marzo 2001 204

Un esempio d’autore

Il seguente problema e dovuto a Hugo Stei-

nhaus puo essere enunciato come segue:

Quattro cani si trovano ai quattro

vertici A,B,C,D di un prato quadra-

to. I cani cominciano a correre con

velocita costante puntando ciascuno

il cane che si trova nel vertice suc-

cessivo. Dopo quanto tempo i quat-

tro cani si incontrano e quanta stra-

da hanno percorso? Quale percorso

hanno compiuto?

29 marzo 2001 205

Per simmetria i quattro cani si trovano sem-

pre ai vertici di un quadrato il cui lato dimi-

nuisce con la velocita con cui ciascun cane si

muove

I quadrati ai vertici dei quali si trovano i cani

in ciascun istante sono centrati nel centro O

di tutti i quadrati

I cani si muovono seguendo una traiettoria

che si mantiene tangente ad uno dei lati dei

quadrati

29 marzo 2001 206

ogni cane punta il successivo cioe il vetto-

re tangente alla curva percorsa (vettore ve-

locita) punta verso il successivo vertice del

quadrato e quindi giace su un lato del qua-

drato.

Ciascun cane si muove seguendo una traiet-

toria per la quale in ogni punto il vettore tan-

gente ed il raggio vettore rispetto al centro

O formano un angolo costante ed uguale aπ4.

Ogni cane raggiunge il successivo nel tempo

necessario a percorrere il lato del quadrato e

quindi percorre una distanza pari al lato del

quadrato stesso.

La traiettoria percorsa e una spirale logarit-

mica centrata in O e relativa ad un angolo

tra tangente e raggio vettore di π4

29 marzo 2001 207

Consideriamo un riferimento polare (ρ, θ) cen-

trato in O e una parametrizzazione della cur-

va descritta dal cane che parte dal vertice in

A ed insegue il cane che inizialmente e in D

In un successivo istante i cani si trovano in

A(t) = A′, B(t) = B′, C(t) = C′, D(t) = D′

e la traiettoria percorsa e caratterizzata dal

fatto che il raggio vettore A′O forma un an-

golo di π4 con il lato A′D′ che risulta essere

tangente alla traiettoria.

29 marzo 2001 208

Siano (x(t), y(t)) le equazioni parametriche

della curva descritta dal cane in A e sia ρ(t) =

ρ(θ(t)); avremo chexa(t) = x0 + ρ(t) cos θ(t)

ya(t) = y0 + ρ(t) sin θ(t)(280)

mentre il cane in D (che e inseguito dal cane

in A) percorrera una traiettoria di equazionixd(t) = x0 + ρ(t) sin θ(t)

yd(t) = y0 − ρ(t) cos θ(t)(281)

Analogamente i cani in B e C percorreranno

rispettivamente traiettorie di equazionixb(t) = x0 − ρ(t) sin θ(t)

yb(t) = y0 + ρ(t) cos θ(t)(282)

xc(t) = x0 − ρ(t) cos θ(t)

yc(t) = y0 − ρ(t) sin θ(t)(283)

29 marzo 2001 209

Il vettore tangente alla curva descritta dalcane in un generico punto A(t) e dato da

T (t) =

{xa(t) = (ρ(t) cos θ(t)− ρ(t) sin θ(t))θ(t)

ya(t) = (ρ(t) sin θ(t) + ρ(t) cos θ(t))θ(t)(284)

mentre il vettore che indica la direzione che

il cane deve seguire e D(t) − A(t) e si puo

scrivere nella forma

D(t)−A(t) =

ρ(t) sin θ(t)− ρ(t) cos θ(t)

−ρ(t) cos θ(t)− ρ(t) sin θ(t)

(285)

Dovra aversi

T (t) = k[D(t)−A(t)] (286)

dove k e una costante di proporzionalita, per

cui

29 marzo 2001 210

(ρ(t) cos θ(t)− ρ(t) sin θ(t))θ(t) =

= k(ρ(t) sin θ(t)− ρ(t) cos θ(t))(ρ(t) sin θ(t) + ρ(t) cos θ(t))θ(t) =

= k(−ρ(t) cos θ(t)− ρ(t) sin θ(t))(287)

Dividendo membro a membro 287 otteniamo

ρ(t) cos θ(t)− ρ(t) sin θ(t)

ρ(t) sin θ(t) + ρ(t) cos θ(t)=

=ρ(t) sin θ(t)− ρ(t) cos θ(t)

−ρ(t) cos θ(t)− ρ(t) sin θ(t)(288)

moltiplicando per i denominatori

− ρρ cos2 θ+ ρ2 sin θ cos θ −− ρρ cos θ sin θ+ ρ2 sin2 θ =

ρρ sin2 θ − ρρ sin θ cos θ+

+ ρ2 sin θ cos θ − ρ2 cos2 θ (289)

29 marzo 2001 211

da cui deduciamo che

ρρ− ρ2 = 0 ρ = ρ (290)

(ovvero

ρρ− ρ2θ(t) = 0 (291)

nel caso in cui

ρ(t) = ρ(θ(t))) (292)

ed infine

ρ(t) = ceθ(t) (293)

La 289 definisce l’equazione polare della tra-

iettoria percorsa dal cane.

Per ottenerne la legge oraria occorre tenere

anche conto del fatto che la curva e percorsa

con velocita costante |v|

29 marzo 2001 212

In altre parole√x2(t) + y2(t) = |v| (294)

da cui

(ρ2(t) + ρ2(t))θ2(t) = |v|2 (295)

Se ne deduce che

2c2e2θ(t)θ2(t) = |v|2 (296)

eθ(t)θ(t) = ±|v|√2c

(297)

e quindi

eθ(t) = ±1√2c|v|t+ b (298)

θ(t) = ln(b±1√2c|v|t) (299)

29 marzo 2001 213

Curva di inseguimento

Se in generale il bersaglio T si muove lungo

una traiettoria γ di equazioni

γ(t) = (a(t), b(t)) (300)

E se P e l’inseguitore che procede a velo-

cita costante v lungo la curva di equazioni

parametriche

(x(t), y(t))

dovra aversix(t) = k(t)(a(t)− x(t))

y(t) = k(t)(b(t)− y(t))

(x(t))2 + (y(t))2 = v2(301)

29 marzo 2001 214

Dalla 301 si ricava che

k(t) =v√

(a(t)− x(t))2 + (b(t)− y(t))2

(302)

e otteniamo per sostituzionex(t) = v a(t)−x(t)√

(a(t)−x(t))2+(b(t)−y(t))2

y(t) = v b(t)−y(t)√(a(t)−x(t))2+(b(t)−y(t))2

(303)

Il sistema puo essere integrato numericamen-

te.

29 marzo 2001 215

Curva di pedinamento

Per il pedinamento, se T percorre γ di equa-

zioni

γ(t) = (a(t), b(t)) (304)

e P e il punto inseguitore, si muove lungo la

curva

(x(t), y(t))

Si avrax(t) = k(t)(a(t)− x(t))

y(t) = k(t)(b(t)− y(t))

(a(t)− x(t))2 + (b(t)− y(t))2 = d2

(305)

29 marzo 2001 216

Dalla 305 si ricava che

2[(a(t)− x(t)][a(t)− x(t)] +

+ [b(t)− y(t)][b(t)− y(t)] = 0 (306)

e quindi

[a(t)− x(t)][a(t)− k(t)(a(t)− x(t))] +

+ [b(t)− y(t)][b(t)− k(b(t)− y(t))] = 0

(307)

e ancora

a(t)[a(t)− x(t)]− k(t)[a(t)− x(t)]2 +

+ b(t)[b(t)− y(t)]− k(t)[b(t)− y(t)]2 = 0

(308)

29 marzo 2001 217

Infine

a(t)(a(t)− x(t) + b(t)− k(t)d2 = 0

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

d2(309)

e

x(t) = a(t)(a(t)−x(t)+b(t)−k(t)d2

(a(t)− x(t))

y(t) = a(t)(a(t)−x(t)+b(t)−k(t)d2

(b(t)− y(t))

(310)

Il sistema puo essere integrato numericamen-

te.

29 marzo 2001 218

La curva di inseguimento nel caso in cui

il bersaglio si muova di moto circolare

29 marzo 2001 219

La curva di inseguimento nel caso in cui il bersaglio

si muova di moto circolare per t = .4 e t = .8

29 marzo 2001 220

La curva di inseguimento

nel caso in cui il bersaglio

si muova di moto circolare

per t = 1.2 e t = 1.6

29 marzo 2001 221

La curva di inseguimento

nel caso in cui il bersaglio

si muova di moto circolare

per t = 2 e t = 2.4

29 marzo 2001 222

La curva di pedinamento

nel caso in cui il bersaglio

si muova di moto rettilineo

29 marzo 2001 223

La curva di pedinamento

nel caso in cui il bersaglio

si muova di moto circolare uniforme

29 marzo 2001 224

La curva di pedinamento

nel caso in cui il bersaglio

si muova di moto circolare

(nel secondo caso la distanza di pedinamento

e uguale al raggio del cerchio

29 marzo 2001 225

La spirale logaritmica

Stabiliamo un riferimento polare nel piano;

O e il polo di tale sistema

ρ e il raggio vettore

θ e l’angolo formato dal raggio vettore con il

semiasse positivo delle x.

29 marzo 2001 226

Un punto e individuato dalle coordinatex = ρ cos θ

y = ρ sin θ(311)

Se e data ρ = ρ(θ) si individua nel piano una

curva di equazioni parametrichex = ρ(θ) cos θ

y = ρ(θ) sin θ(312)

Si definisce spirale logaritmica una curva nel

piano tale che in ogni punto il raggio vettore

della curva ed il vettore tangente alla curva

formino un angolo fissato α

29 marzo 2001 227

Poiche il prodotto scalare dei vettori (x, y) e

x, y) e uguale al prodotto delle norme dei vet-

tori stessi e del coseno dell’angolo compreso.

una curva che soddisfi la proprieta richiesta

deve soddisfare

x(t)x(t) + y(t)y(t) =

= ‖(x(t), y(t))‖‖(x(t), y(t))‖ cosα (313)

Se teniamo conto chex(t) = ρ(θ(t)) cos θ(t)

y(t) = ρ(θ(t)) sin θ(t)(314)

29 marzo 2001 228

e chex(t)x(t) + y(t)y(t) = d

dt12ρ

2(t) = ρ(t)ρ(t)θ(t)

‖(x(t), y(t))‖2 = ρ2(t)

‖(x(t), y(t))‖2 = (ρ2(t) + ρ2(t))θ2(t)(315)

Otteniamo

ρ2(t)ρ2(t)θ2(t) = cos2αρ2(t)θ2(t)(ρ2(t)ρ2(t))

(316)

Pertanto

ρ2(t) = cos2α(ρ2(t) + ρ2(t)) (317)

ρ2(t) =cos2α

1− cos2αρ2(t) (318)

(319)

ed infine

ρ(t) = cotαρ(t) (320)

29 marzo 2001 229

L’ultima equazione descrive la forma della

spirale logaritmica ma non dice nulla riguardo

al modo in cui essa viene percorsa.

Descrive le proprieta geometriche della spira-

le logaritmica ma non le sue proprieta cine-

matiche.

La velocita e descritta dal vettore (x(t), y(t))

la cui norma e

‖(x(t), y(t))‖2 = x2(t) + y2(t) =

(ρ2(t) + ρ2(t))θ2(t) (321)

possiamo imporre che la spirale venga percor-

sa con velocita costante in modulo chiedendo

che

(ρ2(t) + ρ2(t))θ2(t) = |v| (322)

Le equazioni trovate possono essere facil-

mente risolte per trovare ρ(t) e θ(t).

29 marzo 2001 230

Qualche risultato

di

teoria del controllo

Consideriamo un modulo lunare ad altezza

h dal suolo con velocita v0 diretta verso il

basso,

Occorre determinare la minima quantita di

carburante u necessaria per alimentare i mo-

tori che forniscono la forza frenante necessa-

ria per toccare il suolo con velocita nulla.

Si suppone la forza frenante pari a ku (pro-

porzionale al consumo di carburante)

Il modulo ha massa m ed e soggetto alla forza

di gravita lunare mg;

29 marzo 2001 231

Se x(t) e la quota del modulo all’istante t

dovra essere

mx = −mg+ ku(t) (323)

Il carburante u(t) fornito ai motori e limitato

dalla disponibilita totale α per cui

0 ≤ u(t) ≤ α (324)

Il carburante controlla la capacita frenante

del modulo che parte dalle condizioni iniziali

x(0) = h x(0) = v0 (325)

29 marzo 2001 232

A meno di ridenominare le costanti in gioco

possiamo riformulare il problema nei seguenti

termini: x(t) = −g+ ku(t)

0 ≤ u(t) ≤ α

x(0) = h

x(0) = v0

(326)

Possiamo anche supporre k = 1 a meno di

ridenominare anche α, e riscrivere l’equazione

in forma di sistema del primo ordine

x(t) = y(t)

y(t) = −g+ u(t)

0 ≤ u(t) ≤ α

x(0) = h

y(0) = v0

(327)

29 marzo 2001 233

Perche il modulo non si sfracelli al suolo,

dovra essere

x(T ) = y(T ) = 0 (328)

dove T e l’istante in cui il modulo tocca il

suolo; il consumo di carburante puo essere

stimato mediante la∫ T0u(t)dt (329)

Occorre determinare u in modo che sia mini-

mo l’integrale in 329 al variare delle funzioni

u che soddisfano le condizioni 327.

Per risolvere il problema occorre il principio

del massimo di Potriagyn

29 marzo 2001 234

Il problema di controllo

Sia dunque

f0 : Rn × Rk → R

sufficientemente regolare (almeno di classe

C1) e sia

F : Rn × Rk → Rn

una funzione vettoriale sufficientemente re-

golare.

Consideriamo i seguenti, di poco differenti,

problemi di controllo ottimo:

29 marzo 2001 235

Determinare il minimo del funzionale

∫ T0f0(x(t), u(t))dt

al variare di x, funzione di stato del

sistema, e di u, funzione di controllo,

tali che:

x(t) = F (x(t), u(t)) (330)

x(0) = x0 (331)

|u(t)| ≤ α (332)

29 marzo 2001 236

Determinare il minimo del funzionale

∫ T0f0(x(t), u(t))dt

al variare di T , tempo finale, di x,

funzione di stato del sistema, e di u,

funzioni di controllo, tali che:

x(t) = F (x(t), u(t)) (333)

x(0) = x0 (334)

|u(t)| ≤ α (335)

29 marzo 2001 237

Determinare il minimo del funzionale

∫ T0f0(x(t), u(t))dt

al variare di T , tempo finale, di x,

funzione di stato del sistema, e di u,

funzioni di controllo, tali che:

x(t) = F (x(t), u(t)) (336)

x(0) = x0 , x(T ) = x1 (337)

|u(t)| ≤ α (338)

29 marzo 2001 238

Possiamo enunciare i seguenti risultati che

forniscono condizioni necessarie per lo stato

ed il controllo ottimale e sono note come

Principio del massimo di Pontryagin.

Distingueremo tre situazioni

• Il caso in cui T sia fissato,

• Il caso in cui sia richiesta una condizione

finale

x(T ) = x1

• Il caso in cui T non sia fissato

29 marzo 2001 239

Nel caso in cui T sia fissato,

Se (x∗, u∗) e una soluzione del problema di

controllo allora esiste

λ0 ∈ R , λ : R → Rn

tali che, posto

H(x, u) = λ0f0(x, u) + 〈λ, F (x, u)〉 =

= λ0f0(x, u) +n∑i=1

λiFi(x, u) (339)

si abbia

x∗(t) = F (x∗(t), u∗(t)) (340)

x∗(0) = x0 (341)

λ(t) = −∇xH(x∗(t), u∗(t)) (342)

H(x∗(t), u∗(t)) = max|u|≤α

H(x∗(t), u) (343)

λ0 = −1 (344)

λ(T ) = 0 (345)

29 marzo 2001 240

Nel caso in cui sia richiesta una

condizione finale

x(T ) = x1

Valgono le precednti condizioni ma puo acca-

dere che λ0 = 0, ma in tal caso si ha λ(t) 6= 0

per ogni t ∈ [0, T ].

Infine

Nel caso in cui T non sia fissato

Oltre alle precedenti condizioni deve essere

verificata la condizione

H(λ(t), x∗(t), u∗(t)) = 0

29 marzo 2001 241

Qualche esempio

di

problema di controllo

29 marzo 2001 242

Il problema

del

minimo tempo

Consideriamo una particella, che

si muove su una retta,

x(t) e la distanza della particella

da un punto fisso O, della retta,

detto origine

Il moto della particella e governato

dall’equazione differenziale

x(t) = u(t) |u(t)| ≤ 1 (346)

29 marzo 2001 243

dove u(t) rappresenta una forza ap-

plicata alla particella.

Vogliamo trovare u(t) in modo chela particella che parte con velocitae posizione iniziale nota, raggiun-ga l’origine con velocita nulla nelminimo tempo possibile.

29 marzo 2001 244

Possiamo porre il problema in termini di con-trollo come segue

Minimizzare∫ T0dt (347)

Sugli archi che soddisfano l’equazione di con-

trollox(t) = y(t)

y(t) = u(t)

x(0) = x0y(0) = v0

(348)

ed i vincoli x(T ) = 0

y(T ) = 0(349)

29 marzo 2001 245

Con riferimento all’enunciato del principio del

massimo avremo che

f0((x, y), u) = 1 F ((x, y), u) =

yu

(350)

per cui

H(x, u) = λ0 + λy+ µu (351)

Per il principio del massimo dovra essere(λµ

)= −

(0λ

)=

(HxHy

)(352)

da cui si ricava cheλ(t) = h

µ(t) = −ht+ k(353)

29 marzo 2001 246

Poiche H e lineare in u, assumera il suo mas-simo per u = ±1

Il segno di u deve essere scelto in dipendenzadel segno di µ.

Piu precisamente

u = 1 nel caso in cui µ > 0

u = −1 nel caso in cui µ < 0

Dal momento che µ(t) = −ht+ k si annullauna ed una sola volta, in t, anche u potraavere un solo cambiamento di segno

Pertanto assumera valore costante

±1 per t ∈ [0, t]

e

∓1 per t ∈ [t, T ]

29 marzo 2001 247

Quando u = 1 avremo chex(t) = y(t)

y(t) = 1(354)

e x(t) = t2/2 + at+ b = (t+a)2

2 + b− a2

2y(t) = t+ a

(355)

Le traiettorie per u = 1 sono come in figura

29 marzo 2001 248

mentre per u = −1 avremo chex(t) = y(t)

y(t) = −1(356)

ex(t) = −t2/2 + at+ b = −(t−a)22 + b+ a2

2y(t) = −t+ a = −(t− a)

(357)

Le traiettorie per u = −1 sono come in figura

29 marzo 2001 249

In ogni caso si avra

x = ±y2 + (b∓a2

2) (358)

Le traiettorie del sistema per u = ±1

nel piano (x, y), x(t) individua posizione e y(t)

la velocita della particella.

29 marzo 2001 250

Per raggiungere l’origine, occorrera allora per-

correre una delle due traiettorie che passano

per (0,0) avendo cura di scegliere il verso

diretto all’origine e non quello dall’origine.

Le traiettorie del sistema

29 marzo 2001 251

A questo proposito osserviamo che la traiet-

toria γ che passa per l’origine e corrisponde

a x positive e percorsa per u = 1 concor-

demente al verso positivo di y in quanto su

di essa y = t+ a mentre la traiettoria δ che

passa per l’origine e corrisponde a x negati-

ve e percorsa per u = −1 in senso opposto

al verso positivo di y in quanto su di essa

y = −t+ a.

Queste osservazioni ci permettono di con-

cludere che per raggiungere l’origine ab-

biamo a disposizione soltanto uno dei due

archi di traiettoria che abbiamo descritto

e che sono indicati nella figura seguente.

29 marzo 2001 252

Le traiettorie che conducono all’origine sono

illustrate nella figura seguente

29 marzo 2001 253

Se le condizioni di partenza (x0, v0) = (x0, y0)

giacciono su uno di tali archi di curva possia-

mo allora trovare una soluzione del proble-

ma, semplicemente tenendo fisso u = ±1, a

seconda dei casi.

In caso contrario dovremo, partendo dal pun-

to iniziale assegnato e percorrendo una delle

altre traiettorie,

x = ±y2 + (b∓a2

2) (359)

raggiungere la parte che conduce all’origine

di γ o di δ.

29 marzo 2001 254

Le traiettorie ottimali sono riportate nella

figura seguente

Cio puo essere facilmente fatto imponendo

che

x0 = ±y20 + (b∓a2

2) (360)

29 marzo 2001 255

il Problema del modulo lunare

Supponiamo di pilotare una navi-

cella, in discesa verso la superficie

della luna da una quota iniziale h

con velocita iniziale v0, e quindi

soggetta alla forza di gravita, in

modo da portarla a quota 0 con

velocita nulla.

Possiamo utilizzare per questo sco-

po dei motori in grado di forni-

re un’azione frenante proporziona-

le al consumo di carburante.

Il problema consiste nel soddisfarele condizioni richieste con il mini-mo di carburante possibile.29 marzo 2001 256

Se u(t) e il consumo di carburante,

e g e l’accelerazione di gravita,

a meno di costanti moltiplicative, il moto del

modulo sara governato dal sistemax(t) = y(t)

y(t) = −g+ u(t)

x(0) = h

y(0) = v0(361)

e sara nostro scopo minimizzare la quantita∫ T0u(t)dt (362)

sotto i vincoli

x(T ) = 0

y(T ) = 00 ≤ u(t) ≤ α (363)

29 marzo 2001 257

Con riferimento all’enunciato del principio del

massimo avremo che

f0((x, y), u) = u F ((x, y), u) =

(y

−g+ u

)(364)

per cui

H(x, u) = λ0u+ λy+ µ(−g+ u) =

(λ0 + µ)u+ λy − µg (365)

Per il principio del massimo dovra essere(λµ

)= −

(0λ

)=

(HxHy

)(366)

da cui si ricava cheλ(t) = c

µ(t) = −ct+ d(367)

29 marzo 2001 258

Anche in questo caso H e lineare in u,

Essa assumera il suo massimo per u = 0

oppure u = α

La scelta del controllo u da applicare dipen-

dera dal segno di λ0 + µ.

Piu precisamente si avra

• u = 0 nel caso in cui λ0 + µ < 0

• u = α nel caso in cui µ > 0.

29 marzo 2001 259

Dal momento che λ0(t) + µ(t) = −ct+ d′ siannulla una ed una sola volta, in t,

u potra avere un solo cambiamento di segno

Pertanto assumera valore costante 0 fino a t

e valore α da t in poi (o viceversa).

Quando u = 0 avremo chex(t) = y(t)

y(t) = −g(368)

e x(t) = −gt2

2 + at+ b

y(t) = −gt+ a(369)

mentre per u = α avremo chex(t) = y(t)

y(t) = −g+ α = k(370)

e x(t) = kt2

2 + at+ b

y(t) = kt+ a(371)

29 marzo 2001 260

Se v0 6= 0 non e possibile riuscire a toccare

il suolo con velocita nulla senza frenare (u =

0),

Quindi occorre utilizzare una delle traiettorie

che corrispondono ad u = α

Piu precisamente occorre usare una delle tra-

iettorie che passano per x = 0 ed y = v =

0.

Affinche in un istante t si abbia

x(t) = y(t) = 0

deve risultare

kt+ a = 0 da cui a = −kt

0 =kt2

2+ at+ b =

kt2

2− kt2 + b = −

kt2

2+ b =

(372)

da cui

b =kt2

2

29 marzo 2001 261

Ne deduciamo chey = k(t− t)

x = kt2

2 − ktt+ kt2

2 = k(t−t)22

(373)

ed ancora

x =y2

2k(374)

La 371 definisce una parabola di cui si consi-

dera solo la parte che giace nel quadrante in

cui x, lo spazio, e y, la velocita, sono positivi.

29 marzo 2001 262

Queste osservazioni ci permettono di con-

cludere che per raggiungere l’origine abbia-

mo a disposizione soltanto l’arco di parabo-

la γ definito dalla 371 che giace nel primo

quadrante.

Nella figura che segue sono riportate due tra-

iettorie di caduta libera e la traiettoria che

passa per l’origine e corrisponde al controllo

u = α

29 marzo 2001 263

Se le condizioni di partenza (x0, v0) = (x0, y0)

giacciono sull’arco di curva γ la soluzione del

problema, consiste nell’applicare il controllo

u = α, cioe frenare al massimo per l’intera

discesa.

In caso contrario sara necessario, partendo

dal punto iniziale assegnato, percorrere una

delle traiettorie definite da 369, raggiunge-

re un punto dell’arco γ e da liı proseguire

sull’arco γ.

29 marzo 2001 264

Nel secondo caso il comportamento da te-

nere sara lasciare il modulo in caduta libera

sino a raggiungere un punto della traiettoria

γ e quindi frenare per il restante tempo di

discesa.

La traiettoria ottimale e riportata nella figura

seguente

29 marzo 2001 265

Qualche risultato

di

calcolo delle variazioni

Il problema classico di calcolo delle variazioni

consiste nel minimizzare un funzionale inte-

grale che dipende da una funzione e dalla sua

derivata.

Piu precisamente se

f : R2 → R

e regolare e se C1 e lo spazio delle funzioni

x : [a, b] → R

che siano continue con la loro derivata, pos-

siamo porci il problema di trovare

minx∈C1,x(a)=α,x(b)=β

{∫ baf(x(t), x(t))dt

}(375)

29 marzo 2001 266

La curva di lunghezza minima

per due punti fissati

Il problema consiste nel determinare la curva

y(x), nel piano (x, y),

che passa per due punti

P0 = (a, y0) e P1 = (b, y1)

percorrendo la quale si passa dal punto P0 al

punto P1 percorrendo la minima distanza.

La lunghezza di una curva si calcola mediante

la ∫ ba

√1 + (y′(x))2dx

Il problema si riduce ovviamente a trovare

miny∈C1,y(a)=y0,y(b)=y1

{∫ ba

√1 + (y′(x))2dx

}(376)

29 marzo 2001 267

La brachistocrona

Il problema consiste nel trovare la curva da

seguire per minimizzare il tempo di percor-

renza necessario perche un punto materiale

soggetto alla sola forza di gravita si sposti

da un punto A = (a, h), h > 0 ad un punto

B = (b,0).

Indichiamo con x, y(x) la funzione il cui grafi-

co (x, y(x)) e percorso per passare dal punto

A al punto B,

Il punto materiale descrivera una traiettoria

rispetto alla quale la lunghezza d’arco e data

da

s(x) =∫ xa

√1 + (y′(x))2dx

29 marzo 2001 268

Se consideriamo una partizione dell’intervallo

[a, b]

{a = x0 < x1 < x2 < .... < xn = b}

la lunghezza percorsa relativa all’intervallo

[xk, xk+1] e∫ xk+1

xk

√1 + (y′(x))2dx =√

1 + (y′(ck))2(xk+1 − xk) (377)

dove ck ∈ [xk, xk+1] La velocita di caduta

dipende solo dalla quota per il principio di

conservazione dell’energia e si ha

1

2mv2 = mgy e v =

√2gy(x)

29 marzo 2001 269

Nell’intervallo [xk, xk+1] possiamo considera-

re che

y(x) = y(dk) con dk ∈ [xk, xk+1]

per cui

v =√

2gy(dk)

Il tempo tk necessario per percorrere la parte

di curva relativa all’intervallo [xk, xk+1] sara

tk =

√1 + (y′(ck))

2√2gy(dk)

(xk+1 − xk)

ck, dk ∈ [xk, xk+1] (378)

29 marzo 2001 270

Sommando su k e passando al limite, si ot-

tiene

t =∫ ba

√1 + (y′(x))2√

2gy(x)dx (379)

Il problema e ridotto a minimizzare il fun-

zionale definito nella 379 sulle funzioni y di

classe C1 tali che y(a) = h e y(b) = 0.

miny∈C1,y(a)=h,y(b)=0

∫ ba

√1 + (y′(x))2√

2gy(x)dx

29 marzo 2001 271

Il solido di rotazione

avente minima superficie

Determinare il solido di rotazione che ha mi-

nima superficie.

In un riferimento cartesiano ortogonale (x, y, z)

consideriamo una curva y(x), x ∈ [a, b] tale

che y(a) = α e y(b) = β

e la superficie da essa generata dalla rotazio-

ne attorno all’asse x

Possiamo parametrizzare la superficie nella

formax = x

y = y(x) cos θ

z = y(x) sin θ

x ∈ [a, b] (380)

29 marzo 2001 272

si ha

∂x, y, z

∂x, ρ, θ=

(1 y′(x) cos θ y′(x) sin θ0 −y(x) sin θ y(x) cos θ

)(381)

per cui l’area della superficie si puo calcolare

mediante la∫ bay(x)

√1 + (y′(x))2dx (382)

Il problema si riduce a quello di trovare il

minimo del funzionale

2π∫ bay(x)

√1 + (y′(x))2dx (383)

sotto le condizioni

y(a) = α y(b) = β (384)

miny∈C1,y(a)=α,y(b)=β

∫ bay(x)

√1 + (y′(x))2dx

29 marzo 2001 273

Il solido di rotazione

avente minima resistenza

all’avanzamento

Determinare il solido di rotazione che offre

minima resistenza all’avanzamento.

La resistenza all’avanzamento in un gas di

densita tale da consentire l’ipotesi che le mo-

lecole si riflettano liberamente puo essere cal-

colata considerando la forza generata dagli

urti delle molecole sulla superficie stessa.

Se S e una porzione di superficie piana in-

clinata di un angolo θ rispetto alla direzione

di avanzamento, su S agisce una forza F ta-

le che F∆t deve essere uguale alla differen-

za della quantita di moto delle molecole che

urtano la superficie e rimbalzano via.

F∆t = ∆Q

29 marzo 2001 274

Una particella di massa µ che urta la superfi-

cie S con velocita vi = (v,0,0) uguale alla ve-

locita di avanzamento e si allontana con una

velocita vr uguale in modulo a v, nel piano

individuato da v e dalla normale N alla su-

perficie formando con essa un angolo uguale

all’angolo di incidenza.

Pertanto

vr = (v cos 2θ, v sin 2θ,0)

29 marzo 2001 275

La differenza tra le quantita di moto iniziale

e finale e

∆Q = µvi−µvr = µ(v(1−cos 2θ),−v sin 2θ,0)

Nel tempo ∆t sulla superficie S agisce una

forza F tale che

F∆t = ∆Q

Osserviamo che F ha una componente op-

posta alla direzione del moto ed una compo-

nente ad essa perpendicolare.

Quest’ultima non sortisce effetto in quanto

il solido e di rotazione e la somma di tali

componenti e pertanto nulla.

29 marzo 2001 276

Sulla superficie S agisce quindi complessiva-

mente una forza frenante

F =1

∆tm(v(1− cos 2θ),0,0)

dove m e la massa totale delle particelle di

gas che collidono con la superficie S nel tem-

po ∆t.

Avremo che

m = ρS sin(θ)v∆t

in quanto

29 marzo 2001 277

• nel tempo ∆t collidono con la superficie

S le particelle contenute nel volume di un

cilindro di altezza ` = v∆t e di base S

• se ρ e la densita del gas e V e il volume

considerato

m = ρV

• Il fattore sin θ rende conto dell’inclinazio-

ne di S

29 marzo 2001 278

Pertanto su una superficie S inclinata di un

angolo θ agisce una forza pari a

1

∆tv(1−cos(2θ))ρS sin(θ)v∆t = 2ρv2S sin3(θ)

Supponiamo che il solido di rotazione sia de-

finito in un riferimento cartesiano ortogonale

(x, y, z) mediante la rotazione attorno all’asse

x della curvax = x x ∈ [a, b]

y = y(x)(385)

con y(0) = 0 e y(a) = R.

29 marzo 2001 279

Una parametrizzazione s di S e data da

s(x, t) =

x = x

y = y(x) cos t

z = y(x) sin t

x ∈[0, a]

t ∈[0,2π](386)

Si ha

∇s(x, t) =

(1 y′(x) cos t y′(x) sin t0 −y(x) sin t y(x) cos t

)(387)

e

N = (y(x)y′(x),−y(x) cos t,−y(x) sin t)

(388)

da cui

dσ = y(x)√

(1 + y′(x))2 (389)

29 marzo 2001 280

La forza agente su un elemento ∆Sdi super-

ficie sara

2ρv2∆S sin3(θ)

per cui la resistenza si puo calcolare mediante

la∫S2ρv2 sin3(θ)dσ =∫ a

0

∫ 2π

02ρv2 sin3(θ)y(x)

√1 + (y′(x))2dxdt =

4πρv2∫ a0

sin3(θ)y(x)√

1 + (y′(x))2dx (390)

Poiche

sin θ =y′(x)√

1 + (y′(x))2

La forza resistente sara

4πρv2∫ a0

(y′(x))3y(x)

1 + (y′(x))2dx

29 marzo 2001 281

Il problema consiste nel trovare il minimo del

funzionale

4πρv2∫ a0

(y′(x))3y(x)

1 + (y′(x))2dx

sulle funzioni y ∈ C1 tali che , y(0) = 0 e

y(a) = R cioe trovare

miny∈C1,y(0)=0,y(a)=R

4πρv2∫ a0

(y′(x))3y(x)

1 + (y′(x))2dx

29 marzo 2001 282

Le condizioni necessarie

per

il minimo di un funzionale

Sia

f : R2 → R

continua con la sua derivata

Consideriamo lo spazio C1 delle funzioni

x : [a, b] → R

continue con la loro derivata;

ed il problema di trovare

minx∈C1,x(a)=α,x(b)=β

∫ baf(x(t), x(t))dt (391)

29 marzo 2001 283

Supponiamo cioe che

x0 ∈ C1 , x0(a) = α , x0(b) = β

sia tale che∫ baf(x0(t), x0(t))dt =

minx∈C1,x(a)=α,x(b)=β

∫ baf(x(t), x(t))dt (392)

Sia h ∈ C1 tale che h(a) = h(b) = 0

la funzione

φ(λ) =∫ baf(x0(t) + λh(t), x0(t) + λh(t))dt

ammette minimo in λ = 0 per cui dovra

aversi

φ′(0) = 0

29 marzo 2001 284

D’altro canto

φ′(0) =

limλ→0

1

λ

∫ b

a

(f(x0(t) + λh(t), x0(t) + λh(t))− f(x0(t), x0(t)

)dt =

=

∫ b

a

(fx(x0(t), x0(t))h(t) + fv(x0(t), x0(t))h(t)

)dt =

=

∫ b

a

fx(x0(t), x0(t))h(t)dt+ fv(x0(t), x0(t))h(t)∣∣∣ba−

−∫ b

a

(d

dtfv(x0(t), x0(t))

)h(t)dt =

=

∫ b

a

(fx(x0(t), x0(t))−

d

dtfv(x0(t), x0(t))

)h(t))dt = 0

(393)

e si puo dedurre che deve essere

fx(x0(t), x0(t))−d

dtfv(x0(t), x0(t)) = 0

(394)

29 marzo 2001 285

La 394 e l’ equazione di Eulero o euleriana

del problema

Se l’integranda f e autonoma, la 394 puo

essere ulteriormente elaborata come segue

fx(x(t), x(t))−d

dtfv(x(t), x(t)) =

fx(x(t), x(t))− fv,x(x(t), x(t))x(t)−− fv,v(x(t), x(t))x(t) (395)

e moltiplicando per x(t)

fx(x(t), x(t))x(t)− fv,x(x(t), x(t))x2(t)−

− fv,v(x(t), x(t))x(t)x(t) =

= fx(x(t), x(t))x(t)+fv(x(t), x(t))x(t)−fv(x(t), x(t))x(t)−− fv,x(x(t), x(t))x

2(t)− fv,v(x(t), x(t))x(t)x(t) =

=d

dt

(f(x(t), x(t))− x(t)fv(x(t), x(t))

)= 0 (396)

29 marzo 2001 286

Ne deduciamo che perche la 394 sia soddi-

sfatta e necessario e sufficiente che valga

f(x(t), x(t))− x(t)fv(x(t), x(t)) = costante

(397)

L’ultima equazione si chiama integrale pri-

mo della euleriana e spesso e piu facilmente

risolubile della euleriana originale.

29 marzo 2001 287

La soluzione dei problemi

29 marzo 2001 288

La curva di lunghezza minima

per due punti fissati

Avremo

f(y, v) =√

1 + v2

e pertanto

fy(y, v) = 0 fv(y, v) =v√

1 + v2

e la 394 puo essere scritta nella forma

d

dx

y′(x)√1 + (y′(x))2

= 0 (398)

29 marzo 2001 289

Ne viene

y′′(x)√

1 + (y′(x))2 − y′(x) y′(x)y′′(x)√1+(y′(t))2

1 + (y′(x))2=

= y′′(x)1

(1 + (y′(x))2)√

1 + (y′(x))2= 0

(399)

e pertanto se ne ricava

y′′(x) = 0

e

y(x) = ax+ b

Le costanti a, b si determinano imponendo il

passaggio per i punti assegnati.

29 marzo 2001 290

La brachistocrona

Avremo

f(y, v) =

√1 + v2√y

e dalle 397 possiamo ricavare che deve essere√1 + (y′(x))2√

y(x)− y′(x)

y′(x)√y(x)

√1 + (y′(x))2

= c

(400)

Se ne ricava che

1√y(x)

√1 + (y′(x))2

= c (401)

ed elevando al quadrato,

y(x) =k

(1 + (y′(x))2)(402)

29 marzo 2001 291

Procedendo come indicato nel paragrafo de-

dicato alle equazioni differenziali di tipo par-

ticolare e tenendo conto che nel nostro caso

avremo che

f(p) =k

1 + p2

avremo una soluzione in corrispondenza delle

seguenti equazioni parametrichex(t) = x0 + k

∫ tt0

−2(1+φ2(s))2

ds

y(t) = k1+φ2(t)

(403)

Si ottiene una notevole semplificazione dei

calcoli se si sceglie

φ(t) = tan t

in tal caso infattix(t) = x0 − k

∫ tt0

2cos2 sds = x0 − k∫ tt0(1− cos 2s)ds =

= x0 − k(t− 2 sin 2t)

y(t) = cos2 t = 1−cos(2t)2

(404)

29 marzo 2001 292

Il solido di rotazioneavente minima superficie

Avremo

f(y, v) = y√

1 + v2

e dalle 397 possiamo ricavare che deve essere

y(x)√

1 + (y′(x))2 − y′(x)y(x)y′(x)√1 + (y′(x))2

= c

(405)

da cui

y(x)(1 + (y′(x))2 − (y′(x))2y(x)√1 + (y′(x))2

= c (406)

e

y(x)√1 + (y′(x))2

= c (407)

da cui

y(x) = c√

1 + (y′(x))2 (408)

29 marzo 2001 293

Procedendo come indicato nel paragrafo de-

dicato alle equazioni differenziali di tipo par-

ticolare e tenendo conto che nel nostro caso

avremo che

f(p) = c√

1 + p2

avremo una soluzione in corrispondenza delle

seguenti equazioni parametrichex(t) =

∫ tt0

cφ′(s)√1+φ2(s)

ds

y(t) = c√

1 + φ2(t)(409)

ed e evidente che si ottiene una notevole

semplificazione dei calcoli se si sceglie

φ(t) = sinh t

in tal caso infattix(t) =∫ tt0c cosh(s)cosh(s) ds =

∫ tt0cds = ct− ct0

y(t) = sinh(t)(410)

29 marzo 2001 294

Il solido di rotazione

avente minima

resistenza all’avanzamento

Avremo

f(y, v) =v3y

1 + v2

e dalle 397 possiamo ricavare che deve essere

(y′(x))3y(x)

1 + (y′(x))2− y(x)

(y′(x))3(3 + (y′(x))2)

(1 + (y′(x))2)2=

=(y′(x))3y(x)

(1 + (y′(x))2)2(−2) = Costante (411)

Se ne ricava che deve essere

(y′(x))3y(x)

(1 + (y′(x))2)2(−2) = C (412)

ed anche

y(x) = C(1 + (y′(x))2)2

(y′(x))3(413)

29 marzo 2001 295

Possiamo pertanto procedere come indicato

nel paragrafo dedicato alle equazioni diffe-

renziali di tipo particolare ed ottenere una

soluzione in forma parametrica.

In questo caso si ha

f(p) =(1 + p2)2

p3

e possiamo concludere chex(t) = C∫ tt0

(1s −

2s3− 3s5

)ds

y(t) = C (1+s2)2

s3

(414)

Poiche si verifica subito che y(t) non puo mai

assumere il valore nullo, si vede che il pro-

blema nella sua formulazione completa non

ammette soluzione.

29 marzo 2001 296

Procediamo ad alcune approssimazioni

Nel caso in cui y′ sia molto grande, possiamo

supporre che

f(y, v) =v3y

1 + v2= vy

v2

1 + v2≈ vy

ed in tal caso∫ a0

(y′(x))3y(x)

1 + (y′(x))2dx ≈

∫ a0y′(x)y(x)dx =

=(y′(a))2

2= R2 (415)

e quindi il funzionale e indipendente dalla

scelta della funzione.

Nel caso in cui y′ sia trascurabile possiamo

supporre che

f(y, v) =v3y

1 + v2= y(v3 − v5 + v7 + ...) ≈ v3y

29 marzo 2001 297

ed in tal caso∫ a0

(y′(x))3y(x)

1 + (y′(x))2dx ≈

∫ a0

(y′(x))3y(x)dx

Dalle 397 possiamo ricavare che deve essere

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

= −2(y′(x))3y(x) = C (416)

Ne viene che

y′(x)) =C

3√y(x)

(417)

la cui soluzione e data da

y(x) = (αx+ β)34 (418)

Imponendo le condizioni al contorno si ricava

che

y(x) = R

(x

a

)34

29 marzo 2001 298

Qualche risultato

sulle

equazioni differenziali

29 marzo 2001 299

Equazioni differenziali

a variabili separabili

e

Lemma di Gronwall

29 marzo 2001 300

Equazioni differenziali a variabili separabili

Risolvere una equazione differenziale a varia-

bili separabili, significa trovare una funzione

y, derivabile per cui si abbia

y′(x) = f(x)g(y(x)) (419)

con f, g assegnate.

29 marzo 2001 301

Piu precisamente per

I, J ⊂ R

intervalli aperti e non vuoti ed

f : I → R , g : J → R

funzioni continue,

diciamo che y(x) risolve l’equazione differen-

ziale a variabili separabili

y′(x) = f(x)g(y(x)) (420)

se la 420 e soddisfatta da y(x) su un inter-

vallo aperto I ′ ⊂ I.

29 marzo 2001 302

Se y e tale che

g(y) = 0

y(x) = y

e soluzione costante

se g 6≡ 0 per la continuita possiamo sempre

supporre

g(x) 6= 0 ∀x ∈ J (421)

e riscrivere la 420 come

y′(x)

g(y(x))= f(x) (422)

Se F ′ = f in I J e G′ = 1/g in J otteniamo

29 marzo 2001 303

G(y(x)) = F (x) + c , c ∈ R (423)

Ora G e invertibile in J in quanto G′ = 1/g 6=0 ed inoltre si puo sempre scegliere c in modo

che

R(G) ∩R(F + c) = D(G−1) ∩R(F + c) 6= ∅(424)

(cio e necessario affinche la 423 abbia solu-

zione)

Infatti se x0 ∈ I ed y0 ∈ J e c = G(y0)−F (x0),

per la continuita di F + c e G,

R(F + c) e R(G) sono intervalli aperti conte-

nenti il punto G(y0) = F (x0) + c;

Ne concludiamo R(G)∩R(F + c) = D(G−1)∩R(F + c) e a sua volta un intervallo aperto

non vuoto (contiene G(y0)).

29 marzo 2001 304

Se I ′ e l’intervallo degli x tali che la

G(y(x)) = F (x) + c , c ∈ R (425)

e vera

y(x) = G−1(F (x) + c) ∀x ∈ I ′ (426)

Osservazione due primitive F e G possono

essere ottenute nella forma.

G(y) =∫ yy0

ds

g(s)F (x) =

∫ xx0f(t)dt

(427)

nel qual caso si ha

G(y0) = F (x0) = c = 0

29 marzo 2001 305

Il Lemma di Gronwall

Il lemma di Gronwall risolve una disequazio-

ne differenziale che consente stime a priori

per le soluzioni di equazioni differenziali.

Siano

y, f : I → R

continue e positive su un intervallo I, tali che

y(x) ≤∣∣∣∣∣∫ xx0f(t)y(t)dt

∣∣∣∣∣+ c ∀x ∈ I, c > 0

(428)

il Lemma di Gronwall asserisce che

0 ≤ y(x) ≤ ce|∫ xxof(t)dt| ∀x ∈ I (429)

Supponiamo per semplicita x ≥ x0 ;

29 marzo 2001 306

dividendo ambo i membri per il secondo e

moltiplicando poi per f(x) si ottiene

y(x)f(x)

c+∫ xx0f(t)y(t)dt

≤ f(x) (430)

da cui

d

dx

[ln

(c+

∫ xx0f(t)y(t)dt

)]≤ f(x). (431)

integrando

ln

(c+

∫ xx0f(t)y(t)dt

)− ln c ≤

∫ xx0f(t)dt

(432)

onde

c+∫ xx0f(t)y(t)dt ≤ ce

∫ xxof(t)dt (433)

e dalla 428

y(x) ≤ c+∫ xx0f(t)y(t)dt ≤ ce

∫ xxof(t)dt (434)

29 marzo 2001 307

Osserviamo anche che se

y(x) ≤∣∣∣∣∣∫ xx0f(t)y(t)dt

∣∣∣∣∣ ∀x ∈ I (435)

si ha

y(x) ≤∣∣∣∣∣∫ xx0f(t)y(t)dt

∣∣∣∣∣+ c ∀c > 0 (436)

e pertanto

0 ≤ y(x) ≤ ce|∫ xxof(t)dt| ∀c > 0

(437)

per cui, al limite per c→ 0+, si ha y(x) ≡ 0.

29 marzo 2001 308

Le equazioni autonome del primo ordine

Si dicono autonome le equazioni il cui secon-

do membro non dipende esplicitamente dalla

variabile indipendente.

y′(x) = f(y(x)) (438)

Spesso la variabile indipendente rappresenta

il tempo

In tal caso una equazione autonoma model-

lizza fenomeni il cui svolgimento e indipen-

dente dal momento in cui si svolgono.

(Un pendolo oscillera nello stesso modo oggi

come domani e altrettanto indipendente dal

momento del lancio e il moto di una pietra

lasciata cadere da una certa altezza).

29 marzo 2001 309

Per le equazioni autonome le traslate di una

soluzione sono esse stesse soluzioni.

Infatti se y(x) risolve la 438 allora, definito

z(x) = y(x− α), avremo che

z′(x) = y′(x− α) = f(y(x− α)) = f(z(x))

(439)

e quindi z risolve 438.

29 marzo 2001 310

Se f si annulla cioe se f(y) = 0,

y(x) = y e una soluzione (costante dell’equa-

zione 438.

Le soluzioni costanti sono, nel caso delle equa-

zioni autonome, particolarmente significati-

ve.

Possiamo formulare un problema di Cauchy

autonomo come seguey′(x) = f(y(x))

y(x0) = y0(440)

e possiamo osservare che, poiche le traslate

di una soluzione sono esse stesse soluzioni,

non e restrittivo supporre che x0 = 0.

29 marzo 2001 311

E facile risolvere l’equazione autonoma 438 o

equivalentemente il problema di Cauchy 440

al variare dei dati iniziali x0 ed y0, o me-

glio per quanto osservato precedentemente

al variare di y0 e per x0 = 0 (le soluzioni

trascurate si ottengono semplicemente per

traslazione nel tempo delle soluzioni trovate).

Si ha infatti ∫ y(x)y0

1

f(s)ds = x (441)

e si possono ottenere informazioni su y(x)

studiando la sua inversa x(y) definita dalla

∫ yy0

1

f(s)ds = x(y) (442)

29 marzo 2001 312

Osservazioni

Tranne che nelle vicinanze delle soluzioni co-

stanti, f(y) 6= 0 e quindi 1f ha segno costante

ed il primo membro della 442 e certamente

invertibile almeno localmente.

Quando y = y l’integrale che figura a pri-

mo membro della 442 deve essere considera-

to in senso improprio e quindi sara necessario

studiarne la convergenza.

Se l’integrale in questione e convergente la

funzione integrale assume valore x finito in y

e la sua derivata 1f(y) tende a +∞ o a −∞.

Percio in x la soluzione puo essere ”attac-

cata” alla soluzione costante y(x) = y e si

possono (non pero si debbono) verificare fe-

nomeni di non unicita della soluzione.

29 marzo 2001 313

Se invece l’integrale e divergente la soluzione

e definita almeno su una semiretta (non limi-

tata a destra o a sinistra) e risulta asintotica

alla soluzione costante.

Se ad esempio f ∈ C1 la soluzione esiste ed

e unica, inoltre l’integrale in questione e di-

vergente e ogni soluzione e asintotica ad una

soluzione costante.

Gli eventuali zeri della funzione f (soluzioni

costanti) dividono il piano (x, y) in strisce in

cui f e di conseguenza y′ ha segno costante.

Cio permette di ottenere un ragionevole stu-

dio qualitativo delle soluzioni dell’equazione

in modo semplice

29 marzo 2001 314

Due equazioni in forma non normale

Studiamo una coppia di equazioni differen-

ziali in forma non normale che compaiono

nelle condizioni necessarie per la soluzione di

qualche problema di calcolo delle variazioni.

29 marzo 2001 315

Equazioni della forma y(x) = f(y′(x))

Consideriamo una equazione nella forma

y(x) = f(y′(x)) (443)

ove f e una funzione sufficientemente rego-

lare (ad esempio di classe C1. Cerchiamo so-

luzioni derivabili due volte e deriviamo ambo

i membri dell’equazione. Avremo

y′(x) = f ′(y′(x))y′′(x) (444)

cioe

1 =f ′(y′(x))y′′(x)

y′(x)(445)

da cui, integrando,

x = x0 +∫ xx0

f ′(y′(t))y′′(t)

y′(t)dt =

∫ y′(x)y′(x0)

f ′(s)

sds

(446)

Posto y′(x) = p possiamo allora concludere

che x = x0 +∫ pp0

f ′(s)s ds

y = f(p)(447)

29 marzo 2001 316

Alternativamente possiamo cercare soluzioni

parametriche dell’equazione nella formax = x(t)

y = y(t)(448)

Se poniamo

φ(t) = y′(x(t)) (449)

avremo dalla 443

y(t) = y(x(t)) = f(y′(x(t))) = f(φ(t))

ed anche, poiche y(t) = y(x(t))

y(t) = y′(x(t))x(t)

da cui

x(t) =y(t)

y′(x(t))=f ′(φ(t))

φ(t)φ′(t)

Ne viene che la soluzione del problema puo

essere individuata parametricamente dallex(t) =∫ tt0f ′(φ(s))φ(s) φ′(s)ds

y(t) = f(φ(t))(450)

29 marzo 2001 317

Equazioni della forma x = f(y′(x))

Consideriamo una equazione nella forma

x = f(y′(x)) (451)

dove f e una funzione sufficientemente rego-

lare (ad esempio di classe C1).

Cerchiamo soluzioni derivabili due volte e de-

riviamo ambo i membri dell’equazione. Avre-

mo

1 = f ′(y′(x))y′′(x) (452)

da cui,moltiplicando per y′(x)

y′(x) = y′(x)f ′(y′(x))y′′(x) (453)

ed integrando,

y(x)− y0 =∫ y′(x)y′(x0)

f ′(s)sds (454)

29 marzo 2001 318

Posto y′(x) = p possiamo allora concludere

che

x = f(p)

y = y0 +∫ pp0sf ′(s)ds

(455)

Alternativamente cerchiamo soluzioni para-

metriche della 451 nella forma

x = x(t)

y = y(t)(456)

Posto

φ(t) = y′(x(t)) (457)

29 marzo 2001 319

Avremo dalla 451

x(t) = f(y′(x(t))) = f(φ(t))

e dalla 456

y(t) = y(x(t)) e y(t) = y′(x(t))x(t)

da cui

y(t) = x(t)y′(x(t)) = f ′(φ(t))φ(t)φ′(t)

La soluzione del problema puo essere indivi-

duata parametricamente dallex(t) = f(φ(t))

y(t) =∫ tt0f ′(φ(s))φ(s)φ′(s)ds

(458)

29 marzo 2001 320

Stabilita dei sistemi lineari

29 marzo 2001 321

Richiami sui sistemi

differenziali lineari

Raccogliamo qui i principali risultati della teo-

ria dei sistemi di equazioni differenziali lineari.

29 marzo 2001 322

Sia I ⊂ R intervallo aperto non vuoto, non

necessariamente limitato,

Sia

A : I →Mn

una funzione a valori matrici n× n

B : I → Rn

A = (aij) B = (bj) bj, aij : I → R

supponiamo che le funzioni ai,j e bj siano

continue.

Risolvere il sistema differenziale lineare del

primo ordine

Y ′(x) = A(x)Y (x) +B(x) (459)

significa trovare Y : I → Rn derivabile, tale

che la 459 e verificata.

29 marzo 2001 323

L’insieme T di tutte le soluzioni di 459 si

chiama integrale generale del sistema

Il sistema puo essere riscritto usando le com-

ponenti di Y , A, B, nella seguente maniera

y′1(x)y′2(x)...y′n(x)

=

a11(x) a12(x) . . . a1n(x)a21(x) a22(x) . . . a2n(x)

... ... . . . ...an1(x) an2(x) . . . ann(x)

y1(x)y2(x)

...yn(x)

+

b1(x)b2(x)

...bn(x)

(460)

ed anche, in forma piu compatta,

y′i(x) =n∑

j=1

aij(x)yj(x) + bi(x) , i = 1, ..., n.

29 marzo 2001 324

Se B ≡ 0 il sistema si dice omogeneo e assu-

me la forma

Y ′(x) = A(x)Y (x) (461)

Se n = 1 si riduce ad una sola equazione dif-

ferenziale lineare del primo ordine che, posto

A = (a11) = a e B = b1 = b, si scrive nella

forma

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

Nel caso sia assegnato un dato inizialeY ′(x) = A(x)Y (x) +B(x) , ∀x ∈ IY (x0) = Y0

(462)

si parla di problema di Cauchy

29 marzo 2001 325

Nel caso di una equazione, (n = 1), avremoy′(x) = a(x)y(x) + b(x)

y(x0) = y0(463)

E possibile trovare esplicitamente la soluzio-

ne nella forma

y(x) = e∫ xxoa(t)dt

(y0 +

∫ xx0b(t)e−

∫ txoa(s)ds

dt

)(464)

Possiamo riscrivere la 464 nella forma

y(x) = y0e∫ xxoa(t)dt+e

∫ xxoa(t)dt

∫ xx0b(t)e−

∫ txoa(s)ds

dt

29 marzo 2001 326

Tutte le soluzioni dell’equazione omogenea

(b = 0) si ottengono al variare della costante

moltiplicativa y0,

Tutte le soluzioni dell’equazione completa si

ottengono aggiungendo un termine che non

dipende da costanti arbitrarie.

L’insieme delle soluzioni di 463 e uno spa-

zio lineare affine, cioe uno spazio vettoriale

traslato.

Si puo dimostrare che, se i coefficienti del

sistema sono continui allora esiste una ed una

sola soluzione del problema di Cauchy 462

29 marzo 2001 327

Facendo uso di tale risultato possiamo pro-

vare che

L’integrale generale di un sistema differen-

ziale omogeneo del primo ordine e uno spazio

vettoriale avente dimensione uguale al nume-

ro di equazioni del sistema stesso. mentre

l’integrale generale del sistema completo e

uno spazio lineare affine

Piu precisamente se

S = {Y : Y ′ = AY } T = {Y : Y ′ = AY +B}(465)

si ha che

29 marzo 2001 328

S e uno spazio vettoriale di dimensione n

(uguale all’ordine del sistema) mentre T e

uno spazio lineare affine di dimensione n.

Inoltre se Y e una soluzione del sistema com-

pleto, ogni soluzione del sistema completo e

descritta dalla

T = Y + S

La verifica di tali fatti e basata sulla possibi-

lita di stabilire un isomorfismo

Γ : S → Rn

definendo

Γ(Y ) = Y (x0), x0 ∈ I

29 marzo 2001 329

A proposito di Γ si puo ricordare che

1. e banale verificare la linearita di Γ

2. la surgettivita di Γ dipende dal teorema

di esistenza della soluzione del problema

di Cauchy 462

3. l’iniettivita di Γ dipende dal teorema di

unicita della soluzione del problema di Cau-

chy 462

29 marzo 2001 330

Ogni soluzione di un sistema omogeneo di

ordine n puo essere espressa mediante una

combinazione di n soluzioni linearmente indi-

pendenti del sistema stesso.

Siano esse Y1, ..., Yn e sia (yi)j la componente

j-esima della i-esima soluzione.

Possiamo costruire la matrice

G =

(y1)1 (y2)1 . . . (yn)1(y1)2 (y2)2 . . . (yn)2

... ... . . . ...(y1)n (y2)n . . . (yn)n

che indicheremo spesso come

G =(Y1, Y2, ....., Yn

)(466)

considerando gli Yi come vettori colonna

29 marzo 2001 331

G si chiama matrice fondamentale del siste-

ma assegnato.

Ogni soluzione del sistema potra allora essere

scritta nella forma

Y (x) = G(x)C =n∑

j=1

cjYj , C = (cj) ∈ Rn

ovvero, considerando le componenti,

yi(x) =n∑

j=1

(yj)icj.

29 marzo 2001 332

Siano Y1, Y2, ....., Yn n soluzioni del sistema

omogeneo

Y ′(x) = A(x)Y (x)

Chiamiamo wronskiano associato alle n solu-

zioni assegnate il determinante della matrice

(Y1, Y2, ....., Yn).

In altri termini

W (x) =

∣∣∣∣∣∣∣∣∣(y1)1 (y2)1 . . . (yn)1(y1)2 (y2)2 . . . (yn)2

... ... . . . ...(y1)n (y2)n . . . (yn)n

∣∣∣∣∣∣∣∣∣

Siano Y1,Y2,...,Yn n soluzioni del sistema

459. Sono fatti equivalenti:

1. Y1, ..., Yn sono linearmente indipendenti;

2. W (x) 6= 0 ∀x ∈ I ;

3. ∃x0 ∈ I tale che W (x0) 6= 0

29 marzo 2001 333

Se si conosce la soluzione del sistema omo-

geneo

Y ′(x) = A(x)Y (x) , x0 ∈ I

si puo trovare una soluzione del sistema non

omogeneo

Y ′(x) = A(x)Y (x) +B(x)

nella forma

Z(x) = G(x)λ(x)

Dovra aversi

Z′(x) = A(x)Z(x) +B(x)

29 marzo 2001 334

Pertanto, poiche

Z′(x) = G′(x)λ(x) +G(x)λ′(x) ,

deve essere

G′(x)λ(x)+G(x)λ′(x) = A(x)G(x)λ(x)+B(x)

essendo G una matrice fondamentale G einvertibile e

G(x)λ′(x) = B(x) e λ′(x) = G−1(x)B(x).

Se ne deduce che se

λ(x) =∫ xx0

G−1(t)B(t)dt

Z e soluzione.29 marzo 2001 335

Possiamo cosı concludere che se G e unamatrice fondamentale del sistema omogeneo

Y ′(x) = A(x)Y (x)

Allora l’integrale generale del sistema com-pleto

Y ′(x) = A(x)Y (x) +B(x)

e dato da

Y (x) = G(x)

(C +

∫ xx0

G−1(t)B(t)dt

), C ∈ Rn ,

mentre la soluzione del problema di Cauchye data da

Y (x) = G(x)

(G−1(x0)Y0 +

∫ xx0

G−1(t)B(t)dt

)

29 marzo 2001 336

Stabilita di una soluzionedi un sistema autonomo

Nell’ambito dei sistemi lineari, ma piu in ge-nerale dei sistemi di equazioni differenziali eimportante studiare il comportamento del-le soluzioni per x → +∞ (stabilita di unasoluzione).

E necessario che le soluzioni siano almenodefinite sulla semiretta che ha come estremoil punto iniziale

Nel caso dei sistemi autonomi, possiamo sup-porre che il punto iniziale sia l’origine e chela soluzione sia definita su R+ = [0,+∞).

Piu precisamente interessa stabilire se le solu-zioni si avvicinano, per x grande, ad una solu-zione con particolari caratteristiche, e, a me-no di un cambio di variabili, si puo supporreche tale soluzione sia identicamente nulla.

Nel caso dei sistemi lineari omogenei questaipotesi sara automaticamente verificata.

29 marzo 2001 337

Sia dato il sistema

Y ′(x) = F (Y (x)) (467)

dove F : Rn → Rn e lipschitziana di costanteL cioe sia

‖F (y1)− F (y2)‖ ≤ L‖y1 − y2‖ (468)

e F (0) = 0 cosı che la funzione identicamen-te nulla risolve il sistema

In tali condizioni si ha che la soluzione delproblema di Cauchy

Y ′(x) = F (Y (x))

Y (0) = Y0(469)

e definita per x ∈ R+ per ogni y0 ∈ Rn per iteoremi di esistenza ed unicita in grande.

29 marzo 2001 338

Possiamo allora definire

La soluzione nulla e stabile per il sistema se

∀ε > 0 ∃δε > 0 tale che se ‖Y0‖ < δ si ha

‖Y (x)‖ < ε ∀x ≥ 0

Inoltre

La soluzione nulla e asintoticamente stabile

se e stabile e se esiste δ > 0 tale che per

‖Y0‖ < δ ,

limx→+∞

Y (x) = 0.

29 marzo 2001 339

Stabilita dei sistemi lineari 2× 2;

piano delle fasi

Ci servira principalmente conoscere qualche

risultato sulla stabilita dei sistemi lineari a

coefficienti costanti o dei sistemi linearizzati.

Cominciamo ad esaminare la stabilita dei si-

stemi differenziali lineari omogenei a coeffi-

cienti costanti di due equazioni in due inco-

gnite.

29 marzo 2001 340

Consideriamo quindi il sistema

x(t) = ax(t) + by(t)

y(t) = cx(t) + dy(t)a, b, c, d ∈ R

(470)

Possiamo scrivere il sistema anche usando

notazioni vettoriali come

u(t) = Au(t) (471)

con

u(t) =

(x(t)y(t)

)ed A =

(a bc d

)(472)

29 marzo 2001 341

Se u(t) = (x(t), y(t)) e la soluzione di un

sistema, (anche nel caso nonlineare)

Chiamiamo orbita del sistema la curva

descritta parametricamente dalle equazionix = x(t)

y = y(t)(473)

Le orbite (o traiettorie) si rappresentano nel

piano (x, y) (Piano delle Fasi)

Se il sistema e autonomo per ogni punto del

piano delle fasi passa una ed una sola orbita.

Infatti se ci fossero due orbite (x(t), y(t)) e

(ξ(t), η(t)) del sistema tali che

x(t0) = x

y(t0) = y

ξ(τ0) = x

η(τ0) = y(474)

29 marzo 2001 342

allora

ξ1(t) = x(t− τ0 + t0)

η1(t) = y(t− τ0 + t0)(475)

e ancora soluzione del sistema e soddisfa la

condizione

ξ1(τ0) = x

η(τ0) = y(476)

Ne deduciamo che

ξ(t) = ξ1(t) = x(t− τ0 + t0)

η(t) = η1(t) = y(t− τ0 + t0)(477)

e quindi (ξ, η) ed (x, y) percorrono la stessa

orbita a meno di una traslazione nel tempo.

29 marzo 2001 343

Autovalori e Stabilita

Ogni matrice 2 × 2 puo essere trasformata

mediante una matrice di passaggio P , non

singolare in una matrice canonica C (forma

canonica di Jordan) di forma particolarmen-

te semplice, che consente uno studio diretto

della stabilita del sistema differenziale asso-

ciato.

In altre parole e possibile trovare P matrice

non singolare tale che

C = P−1AP

sia una delle matrici che rientrano nel seguen-

te elenco, con riferimento alle caratteristiche

degli autovalori λ1, λ2 della matrice A.

29 marzo 2001 344

1. Nel caso di autovalori reali e distinti λ1,

λ2, (λ1 00 λ2

)

2. Nel caso di autovalori complessi e coniu-

gati λ1 = α+ iβ, λ2 = α− iβ,(α β−β α

)

3. Nel caso di autovalori reali e coincidenti

λ1 = λ2 = λ,(λ 00 λ

)oppure

(λ 0γ λ

)

29 marzo 2001 345

Il sistema dato puo essere trasformato in un

sistema piu semplice la cui matrice dei coeffi-

cienti rientra in uno dei casi appena elencati,

mediante la trasformazione

u = Pv o equivalentemente v = P−1u

Infatti si ottiene

P v = u = Au = aPv (478)

da cui

v = P−1APv = Cv (479)

29 marzo 2001 346

Se siamo in grado di studiare le orbite del si-

stema v = Cv, otteniamo anche informazioni

sul sistema di partenza in quanto la trasfor-

mazione u = Pv e lineare e non degenere e

provoca quindi solo una deformazione delle

orbite.

Diventa importante studiare il comportamen-

to delle orbite dei sistemi associati alle matri-

ci canoniche che abbiamo appena elencato.

Di seguito sono riportati le immagini degli

assi coordinati e del cerchio unitario mediante

la trasformazione

ξ = ax+ by

η = cx+ dyovvero

(ξη

)=

(a bc d

)(xy

)(480)

29 marzo 2001 347

29 marzo 2001 348

29 marzo 2001 349

Caso 1: autovalori reali e distinti

In tal caso il sistema diventa

ξ = λ1ξ

η = λ2η(481)

Se entrambi gli autovalori sono non nulli pos-

siamo allora determinare le soluzioni nella for-

ma

ξ(t) = c1eλ1t

η(t) = c2eλ2t

(482)

Le traiettorie sono date da y = xλ2λ1 e a se-

conda del segno di λ1 e λ2 e del fatto cheλ2λ1

R 1 si possono avere i seguenti casi:

29 marzo 2001 350

Orbite del sistema lineare canonico

con autovalori reali e distinti entrambi

positivi λ2λ1> 1

29 marzo 2001 351

Orbite del sistema lineare canonico

con autovalori reali e distinti entrambi

positivi λ2λ1< 1

29 marzo 2001 352

Orbite del sistema lineare canonico

con autovalori reali e distinti entrambi

negativi λ2λ1> 1

29 marzo 2001 353

Orbite del sistema lineare canonico

con autovalori reali e distinti entrambi

negativi λ2λ1< 1

29 marzo 2001 354

Orbite del sistema lineare canonico

con autovalori reali e distinti uno po-

sitivo ed uno negativo, λ1 < 0 e

λ2 > 0

29 marzo 2001 355

Orbite del sistema lineare canonico

con autovalori reali e distinti uno po-

sitivo ed uno negativo, λ1 > 0 e

λ2 < 0

29 marzo 2001 356

Se uno degli autovalori e nullo allora il siste-

ma diventa

ξ = λ1ξ

η = 0oppure

ξ = 0

η = λ2η

(483)

e le traiettorie sono semirette parallele ad uno

degli assi il cui verso di percorrenza dipende

dal segno dell’autovalore non nullo.

29 marzo 2001 357

Caso 2: autovalori reali e coincidenti

In tal caso il sistema diventa uno dei dueseguenti

ξ = λξ

η = ληoppure

ξ = λξ

η = γξ+ λη

(484)

Cominciamo a supporre che λ 6= 0

Nel primo caso la soluzione del sistema e

ξ(t) = c1eλt

η(t) = c2eλt

(485)

e le orbite sono semirette che passano perl’origine, hanno coefficiente angolare c2

c1e so-

no percorse con verso che si avvicina o siallontana dall’origine a seconda del segno diλ.29 marzo 2001 358

Orbite del sistema lineare canoni-

co con autovalori reali e coincidenti,

positivi

29 marzo 2001 359

Orbite del sistema lineare canoni-

co con autovalori reali e coincidenti,

negativi

29 marzo 2001 360

Nel secondo caso la soluzione del sistema e

ξ(t) = c1eλt

η(t) = (γc1t+ c2)eλt

(486)

da cui si ricava

t =1λ ln

(ξc1

)η = η(t) = (γ1

λ ln(ξc1

)+ c2

c1)ξ

(487)

e si possono disegnare le curve che rappre-

sentano le orbite del sistema.

29 marzo 2001 361

Orbite del sistema lineare canoni-

co con autovalori reali e coincidenti,

negativi

29 marzo 2001 362

Orbite del sistema lineare canoni-

co con autovalori reali e coincidenti,

positivi

29 marzo 2001 363

Se infine λ = 0 otteniamo che, nel primo ca-

so, ogni punto e un punto singolare, mentre

nel secondo caso le soluzioni sono date da

ξ(t) = c1η(t) = γc1t+ c2

(488)

e le orbite sono ancora rette parallele ad un’as-

se.

29 marzo 2001 364

Caso 3: autovalori complessi e coniugati

In tal caso il sistema diventa

ξ = αξ − βη

η = βξ+ αη(489)

e puo essere studiato mediante un ulteriore

cambio di variabili.

Infatti se poniamo

ξ = ρ cos θ

η = ρ sin θ(490)

29 marzo 2001 365

possiamo ricavare che

ρ2(t) = ξ2(t) + η2(t) tan θ(t) =η(t)

ξ(t)(491)

da cui, derivando si ottiene

2ρ(t)ρ(t) = 2ξ(t)ξ(t) + 2η(t)η(t) (492)

e

(493)

(1 + tan2 θ(t))θ(t) =ξ(t)η(t)− η(t)ξ(t)

ξ2(t)(494)

e tenendo conto del sistema

si deduce che ρ = αρ

θ = β(495)

29 marzo 2001 366

da cui si ricava che

ρ(t) = C1eαt

θ(t) = βt+ C2(496)

Le traiettorie sono spirali o, nel caso in cui

α = 0, cerchi.

Il verso di percorrenza e determinato dal se-

gno di β, mentre le spirali che tendono al-

l’origine si distinguono dalle spirali che se ne

allontanano mediante il segno di α (negativo

o positivo rispettivamente)

29 marzo 2001 367

Orbite del sistema lineare canoni-

co con autovalori complessi coniugati

con parte reale positiva

29 marzo 2001 368

Orbite del sistema lineare canoni-

co con autovalori complessi coniugati

con parte reale negativa

29 marzo 2001 369

Orbite del sistema lineare canoni-

co con autovalori complessi coniugati

con parte reale nulla

29 marzo 2001 370

E’ facile provare il seguente criterio di stabi-

lita per i sistemi lineari.

Sia A(t) una funzione a valori matrici n ×n, continua su [t0,+∞) e consideriamo il

sistema lineare

x(t) = A(t)x(t)

Sia G una matrice fondamentale del siste-

ma; allora la soluzione identicamente nulla e

stabile se e solo se

‖G(t)‖ ≤ K , ∀t ≥ t0

inoltre la soluzione nulla e asintoticamente

stabile se e solo se

limt→+∞

‖G(t)‖ = 0

29 marzo 2001 371

Infatti la soluzione x(t) del problema di Cau-

chy

x(t) = A(t)x(t))

x(t0) = x0(497)

si puo esprimere nella forma

x(t) = G(t)G−1(t0)x0 (498)

e

‖x(t)‖ ≤ ‖G(t)‖ ‖G−1(t0)‖ ‖x0‖.

29 marzo 2001 372

Cio permette di concludere che le condizioni

proposte per la stabilita e l’asintotica stabi-

lita della soluzione nulla sono sufficienti.

Viceversa se la soluzione nulla e stabile,possiamo

affermare che per ‖x0‖ < δ si ha ‖x(t)‖ < 1 e

per la 498 se ‖x0‖ < 1 si ha ‖x(t)‖ < 1δ .

Pertanto ogni soluzione che corrisponde al

dato iniziale x0 = ei, dove con (ei) si indi-

ca una ennupla di versori dello spazio Rn,e limitata in norma; ma tali soluzioni sono

le colonne di una matrice fondamentale che

risulta pertanto limitata.

29 marzo 2001 373

Analogamente, se la soluzione nulla e asinto-

ticamente stabile e se ‖x0‖ < δ, si ha

limt→+∞

‖x(t)‖ = 0.

Per la linearita, tutte le soluzioni, e la matrice

G, saranno allora infinitesime per t→ +∞.

Se la matrice dei coefficienti e costante le

soluzioni sono combinazioni lineari di espo-

nenziali, seni e coseni con coefficienti che

possono essere polinomi di t.

Esaminando il comportamento di tali funzioni

si deduce subito che

29 marzo 2001 374

1. se A ha tutti gli autovalori con parte

reale negativa, allora la soluzione nulla

e asintoticamente stabile

2. se A ha qualche autovalore con parte

reale positiva, allora la soluzione nul-

la non e stabile (ne’ asintoticamente

stabile, ovviamente)

3. se A ha tutti gli autovalori con parte

reale negativa o nulla, e gli autovalo-

ri con parte reale nulla sono sempli-

ci, allora la soluzione e stabile ma non

asintoticamente stabile

4. se A ha autovalori con parte reale ne-

gativa o nulla, e gli autovalori con par-

te reale nulla hanno molteplicita supe-

riore 1, allora la soluzione non e asin-

toticamente stabile ma puo essere sia

stabile che instabile a seconda della

forma che assumono le soluzioni.

29 marzo 2001 375

Se il sistema non e lineare, la stabilita e una

questione piu delicata

Si puo, ad esempio, risolvere utilizzando il

metodo delle funzioni di Lyapunov;

Talvolta e possibile studiare la stabilita di

un sistema non lineare utilizzando informa-

zioni sulla stabilita del sistema lineare otte-

nuto usando lo sviluppo di Taylor del primo

ordine della funzione che compare a secondo

membro. (Sistema Linearizzato) .

E interessante conoscere sotto quali condi-

zioni la stabilita del sistema linearizzato e

sufficiente per la stabilita del sistema origi-

nario.

(stabilita in prima approssimazione)

29 marzo 2001 376

A questo proposito possiamo provare che se

x(t) = Ax(t) + f(x(t)) (499)

dove A e una matrice n × n, S = {y ∈ Rn :

‖y‖ < a} ed f : S → R, e continua e tale che

f(0) = 0 e

limx→0

f(x)

‖x‖= 0

se la soluzione nulla e asintoticamente stabile

per il sistema lineare

x(t) = Ax(t) (500)

allora e asintoticamente stabile anche per il

sistema 499 se invece la soluzione nulla non

e stabile per il sistema lineare

x(t) = Ax(t) (501)

allora non e stabile neanche per il sistema

499

29 marzo 2001 377

Dimostrazione. Se ‖x‖ < σ si ha

‖f(x)‖ ≤α

2K‖x‖;

Possiamo supporre che a < σ.

Sia G la matrice fondamentale principale del

sistema linearizzato x = Ax e sia x la soluzio-

ne di 499 definita in un intervallo massimale

[t0, b); posto

z(t) = G(t− t0)x0 +∫ tt0G(t− s)f(x(s))ds

si ha

z′(t) = G′(t− t0)x0+∫ tt0G′(t− s)f(x(s))ds+

+G(0)f(x(t)) =

= A

(G(t− t0)x0 +

∫ tt0G(t− s)f(x(s))ds

)+

+ f(x(t)) =

= Az(t) + f(x(t)) (502)

e

z(t0) = x0.

29 marzo 2001 378

Tenuto conto del fatto che

x(t) = Ax(t) + f(x(t)) , x(t0) = x0

si ha

(z − x)′(t) = A(z − x)(t) , (z − x)(t0) = 0

da cui

z(t) ≡ x(t).

Pertanto

x(t) = G(t− t0)x0 +∫ tt0G(t− s)f(x(s))ds

Poiche gli autovalori di A hanno parte reale

negativa

‖G(t− s)‖ ≤ Ke−α(t−s) , t0 ≤ s ≤ t.

29 marzo 2001 379

Allora

‖x(t)‖ ≤ Ke−α(t−t0)‖x0‖+∫ tt0

α

2e−α(t−s)‖x(s)‖ds

da cui

eα(t−t0)‖x(t)‖ ≤ K‖x0‖+∫ tt0

α

2eα(s−t0)‖x(s)‖ds

e per il lemma di Gronwall

eα(t−t0)‖x(t)‖ ≤ K‖x0‖eα(t−t0)/2

e

‖x(t)‖ ≤ K‖x0‖e−α(t−t0) ≤ K‖x0‖ , ∀t ≥ t0(503)

29 marzo 2001 380

Pertanto si ha, se 0 < ε < σ e ‖x0‖ < ε/K,

‖x(t)‖ < ε , ∀t ∈ [t0, b)

e poiche la soluzione e prolungabile, deve

essere b = +∞.

Passando al limite per t → +∞ nella 503 si

conclude anche l’asintotica stabilita.

Per verificare la seconda affermazione si puo

provare che se G ha qualche autovalore posi-

tivo si puo scegliere una soluzione x in modo

che x(t) → +∞ quando t → +∞, scegliendo

opportunamente i dati iniziali.

29 marzo 2001 381

Qualche osservazione

sul segno degli autovalori nel caso 2× 2

Consideriamo una matrice

A =

(a bc d

)(504)

Gli autovalori sono le soluzioni dell’equazione

det(A− λI) = det

(a− λ bc d− λ

)=

λ2 − (a+ d)λ+ (ad− bc) (505)

Poniamo

T = a+ d (la traccia della matrice)

D = ad− bc (il determinante della matrice)

29 marzo 2001 382

Possiamo riscrivere la 505 come

λ2 − Tλ+D

e gli autovalori sono

λ =T

2±

√T2

4−D (506)

29 marzo 2001 383

Se D > 0

T2

4>T2

4−D

Nel caso in cui T2

4 −D >≥ 0

∣∣∣∣T2∣∣∣∣ >

√T2

4−D

Gli autovalori λ1, λ2 hanno il segno di T

Nel caso in cui T2

4 −D < 0

Gli autovalori sono complessi e coniugati la

loro parte reale <λ1, <λ2 ha lo stesso segno

di T

29 marzo 2001 384

Se D = 0

λ1 = 0 , λ2 = T

Se D < 0

Si ha T2

4 −D > 0

Gli autovalori hanno uno segno positivo e

l’altro negativo.

29 marzo 2001 385

Qualche risultato

sulle

Equazioni differenziali

alle derivate parziali

29 marzo 2001 386

Equazioni lineari in 3 variabili

Consideriamo l’equazione differenziale alle de-

rivate parziali

P (x, y, z)∂u

∂x(x, y, z) +Q(x, y, z)

∂u

∂y(x, y, z) +

+R(x, y, z)∂u

∂z(x, y, z) =

= Pux +Quy +Ruz = 0 (507)

Risolvere la 507 significa trovare una funzio-

ne u(x, y, z) tale che

Pux +Quy +Ruz

Possiamo riscrivere la 507 dividendo per uz( 6=0) per ottenere

Pux

uz+Q

uy

uz+R = 0 (508)

29 marzo 2001 387

Poiche supponiamo uz 6= 0, l’equazione

u(x, y, z) = 0

definisce localmente

z = φ(x, y)

cioe

u(x, u, z) = c ⇐⇒ z = φ(x, y) (509)

e risulta

φx = −ux

uz, φy = −

uy

uz(510)

Quindi u risolve la 507

Pux +Quy +Ruz = 0

se e solo se la corrispondente φ risolve la 508

Pφx +Qφy = R

29 marzo 2001 388

Osserviamo che u(x, y, z) = 0, cosı come z =

φ(x, y), definisce una superficie S in R3 il

cui vettore normale puo essere determina-

to tenendo conto che una parametrizzazione

locale della superficie S e data da

S(x, y) =

x = x

y = y

z = φ(x, y)

(511)

Avremo

∇S(x, y) =

(1 0 φx0 1 φy

)(512)

e quindi

N = (−φx,−φy,1) =(ux

uz,uy

uz,1)

(513)

per cui, a meno di un eventuale fattore di

proporzionalita

N = (φx, φy,−1) oppure N = (ux, uy, uz)

(514)

29 marzo 2001 389

Se definiamo il campo vettoriale

F (x, y, z) =(P (x, y, z), Q(x, y, z), R(x, y, z

)

La 507 oppure, equivalentemente, la 508 puo

essere messa nella forma

〈N,F 〉 = 0 (515)

Pertanto la soluzione u(x, y, z) di 507, o la so-

luzione φ(x, y) della 508, definiscono una su-

perficie S mediante l’equazione u(x, y, z) = 0

oppure z = φ(x, y) la cui normale deve risul-

tare perpendicolare in ogni punto al campo

vettoriale F .

29 marzo 2001 390

Una linea che in ogni punto abbia tangente

parallela al campo F si dice linea di campo.

Poiche la normale alla superficie S soluzione

e perpendicolare al campo F , le linee di forza

devono giacere sulla superficie soluzione S.

Pertanto possiamo descrivere la superficie so-

luzione mediante le linee di forza del campo.

Le linee di forza del campo si possono trovare

nella forma

(x(t), y(t), z(t))

ed in tal caso saranno definite dax(t) = λP (x(t), y(t), z(t))

y(t) = λQ(x(t), y(t), z(t))

z(t) = λR(x(t), y(t), z(t))

(516)

29 marzo 2001 391

Tuttavia possiamo anche esprimere y = y(x)

e z = z(x) per cui

y(t) = y′(x)x(t) , z(t) = z′(x)x(t)

e la 516 puo essere riscritta comey′(x) = Q(x,y(x),z(x))

P (x,y(x),z(x))

z′(x) = R(x,y(x),z(x))P (x,y(x),z(x))

(517)

Se siamo in grado di trovare due funzioni α

e β tali che α(x, y, z) = c1β(x, y, z) = c2

(518)

che forniscano la soluzione di 517 allora avre-

mo

αx + αyy′ + αzz

′ = 0 (519)

cui

αx + αyQ

P+ αz

R

P= 0 (520)

29 marzo 2001 392

Quindi

Rαx +Qαy +Rαz = 0 (521)

In modo simile

Rβx +Qβy +Rβz = 0 (522)

Pertanto α, β sono soluzioni di 507 e si veri-fica subito che

u = Φ(α, β)

e soluzione della stessa 507

Infatti ux = Φααx + Φββx

uy = Φααy + Φββy

uz = Φααz + Φββz

(523)

Da cui

Pux +Quy +Ruz =

= Φα(Rαx+Qαy+Rαz)+Φβ(Rβx+Qβy+Rβz) = 0

(524)

29 marzo 2001 393

Se viceversa u = u(x, y, z) risolve la 507 alloraPux +Quy +Ruz = 0

Rαx +Qαy +Rαz = 0

Rβx +Qβy +Rβz = 0

(525)

Il sistema 525 ammette una soluzione (P,Q,R)

non identicamente nulla e pertanto deve ri-

sultare

det

ux uy uzαx αy αzβx βy βz

= 0 (526)

Quindi si puo dimostrare (Stampacchia Ana-

lisi 2) che esiste una funzione Ψ tale che

u = Ψ(α, β).

29 marzo 2001 394

Equazioni lineari in 2 variabili

Consideriamo l’equazione

a(x, y)ux(x, y) + b(x, y)uy(x, y) = 0 (527)

u(x, y) e soluzione della 527 se se solo se

a(x, y)ux(x, y) + b(x, y)uy(x, y) = 0

Se supponiamo che uy(x, y) 6= 0 avremo che

u(x, y) = u0

mux(x, y(x)) + uy(x, y(x))y

′(x)

m

y′(x) = −ux(x, y(x))

uy(x, y(x))

m

y′(x) = −a(x, y(x))

b(x, y(x))m

g(x, y) = k (528)

29 marzo 2001 395

Quindi u0 e k sono legati da una dipendenza

funzionale; u0 = φ(k) e si ha

u(x, y) = u0 = φ(k) = φ(g(x, y)) (529)

La soluzione dell’equazione differenziale data

puo pertanto essere espressa mediante la so-

luzione g(x, y) = k della ?? e di una funzione

arbitraria φ

Eventuali dati iniziali, ad esempio u(x,0) =

u0(x) possono essere soddisfatti imponendo

condizioni su φ.

Ad esempio

u0(x) = u(x,0) = φ(g(x,0))

da cui si ricava φ.

29 marzo 2001 396

Equazioni Quasi Lineari in 2variabili

Sono equazioni del tipo

P (x, y, φ)∂φ

∂x+Q(x, y, φ)

∂φ

∂y= R(x, y, u)

(530)

Sono equivalenti all’equazione lineare in 3 va-

riabili

P (x, y, z)∂u

∂x(x, y, z) +Q(x, y, z)

∂u

∂y(x, y, z) +

+R(x, y, z)∂u

∂z(x, y, z) = 0 (531)

nell’incognita u(x, y, z) che risulta legata a φ

mediante la

u(x, u, z) = c ⇐⇒ z = φ(x, y) (532)

29 marzo 2001 397

Risulta inoltre

φx = −ux

uz, φy = −

uy

uz(533)

Quindi φ risolve l’equazione quasi lineare

Pφx +Qφy = R

se e solo se la corrispondente u risolve l’e-

quazione lineare

Pux +Quy +Ruz = 0

29 marzo 2001 398

Qualche

Modello Discreto

29 marzo 2001 399

Discreto

Vuol dire

”Fatto di punti isolati”

in contrapposizione con continuo.

E una caratteristica dei Numeri Naturali Nin contrapposizione con la caratteristica di

continuita dei Numeri Reali R.

Possiamo farcene un’idea pensando a come

si rappresentano R ed N su una retta

x x x x x x x x x xh

29 marzo 2001 400

I Conigli di Fibonacci

Nel Liber Abaci e contenuto il seguente pro-

blema:

Quante coppie di conigli verranno

prodotte in un anno a partire da un’u-

nica coppia, se ogni mese ciascuna

coppia genera una nuova coppia che

si riproduce a partire dal suo secondo

mese di vita?

29 marzo 2001 401

La seguente figura illustra l’andamento del

numero di conigli nei primi 5 mesi

29 marzo 2001 402

Consideriamo la successione Fn del numero

di coppie di conigli presenti in ciascun mese.

Avremo

Fn+1 =(numero di coppie presenti al mese

n+1) =(numero di coppie presenti al mese n)

+ (numero di coppie presenti al mese n− 1)

= Fn + Fn−1

In generale si dovra avere

Fn+1 = Fn−1 + Fn

F1 = F2 = 1(534)

Le precedenti uguaglianze definiscono Fn per

ricorrenza.

29 marzo 2001 403

La successione Fn fu chiamata successione diFibonacci da E.E.A Lucas

Possiamo cercare di esplicitare i valori assuntidai suoi termini:

Cerchiamo x ∈ R tale che

Fn = xn

Perche la 534 sia soddisfatta deve essere ve-rificata sia la regola di ricorrenza

Fn+1 = Fn−1 + Fn (535)

che le condizioni iniziali

F1 = F2 = 1 (536)

29 marzo 2001 404

Affinche la 535 sia soddisfatta dovra allora

risultare

xn+1 = xn + xn−1 (537)

e quindi, se escludiamo la soluzione banale

x = 0,

x2 = x+ 1 (538)

da cui

x =1±

√5

2(539)

29 marzo 2001 405

Pertanto

fn =

(1 +

√5

2

)ne gn =

(1−

√5

2

)n(540)

sono due possibili successioni che soddisfano

la regola di ricorrenza 535

Nessuna delle due soddisfa le condizioni ini-

ziali 536

Osserviamo che anche

afn e bgn (541)

soddisfano la 535, infatti

29 marzo 2001 406

afn+1 = afn−1 + afn (542)

e

bgn+1 = bgn−1 + bgn (543)

ed anche

afn + bgn (544)

soddisfa la 535, infatti sommando 542 ed 543

afn+1 + bgn+1 = afn−1 + bgn−1 + afn + bgn

(545)

E tutto questo per ogni a, b ∈ R

29 marzo 2001 407

Possiamo sfruttare la disponibilita di a e b per

sistemare i dati iniziali, infatti avremo che

1 = F1 = af1 + bg1 = a1 +

√5

2+ b

1−√

5

2(546)

1 = F2 = af2 + bg2 =

a

(1 +

√5

2

)2

+ b

(1−

√5

2

)2

(547)

e, risolvendo il sistema lineare si puo ricavare

che

a = −b =1√5

(548)

29 marzo 2001 408

pertanto

Fn =1√5

((1 +

√5

2

)n−(1−

√5

2

)n)(549)

Fn e nota come la successione di Fibonacci

29 marzo 2001 409

Trasmissione di segnali

Supponiamo di voler trasmettere un segnale

su un canale discreto usando due simboli

S1 S2

di uguale durata

d1 = d2 = t secondi

per un tempo T Potremo trasmettere

m =T

tsegnali

cioe

N(T ) = 2m

messaggi diversi. Pertanto la velocita di

trasmissione si potra calcolare mediante il

rapporto

m

T=

log2 2m

T=

log2N(T )

T

29 marzo 2001 410

Se invece si volessero usare due segnali di

diversa durata

d1 6= d2 ad esempio d1 = 1 d2 = 2

potremmo, trasmettere i seguenti messaggi

Unita di tempo Messaggi1 -2 - -

—3 - —

- - -— -

29 marzo 2001 411

E potremmo calcolare il numero di messaggi

mediante le formule di ricorrenza

N(1) = 1

N(2) = 2

N(T ) = N(T − 1) +N(T − 2)

La velocita di trasmissione si esprime median-

te la

log2N(T )

T

che per T grande si stabilizza attorno al va-

lore 0.69 o, piu precisamente, attorno a

√5− 1

2

29 marzo 2001 412

Prezzi di Equilibrio

tra Domanda e Offerta

E molto comune trovare situazioni del gene-re:

• Una azienda agricola decide l’estensione del ter-reno da dedicare alla coltivazione di un certo or-taggio in base al prezzo di tale ortaggio nell’annocorrente.

• Il successivo anno la quantita di ortaggio prodottae eccessiva e cio causa la diminuzione del suoprezzo

• L’azienda diminuisce l’estensione del terreno de-dicato all’ortaggio, e ne diminuisce la produzioneper l’anno seguente

• La minore disponibilita fa aumentare il prezzo del-l’ortaggio e provochera un nuovo aumento dell’e-stensione del terreno coltivato.

• ..............................

• ..............................

29 marzo 2001 413

Per descrivere mediante un modello la situa-

zione definiamo

Sn e l’offerta (supply) di prodotto nell’anno n

Dn e la domanda di prodotto nell’anno n

pn e il prezzo del prodotto nell’anno n

e supponiamo che

• Esiste una relazione tra domanda Dn e

prezzo pn

Dn = f(pn) (550)

ovviamente f e decrescente (la domanda

diminuisce, all’aumentare del prezzo)

• Esiste una relazione tra offerta Sn+1 e

prezzo pn

Sn+1 = g(pn) (551)

29 marzo 2001 414

Ovviamente g e crescente ( l’offerta aumen-

ta, all’aumentare del prezzo)

• Il prezzo di mercato e determinato dalla

domanda e dall’offerta e le contrattazioni

avvengono al prezzo pn+1 in cui domanda

ed offerta sono uguali cioquando si ha

Sn+1 = Dn+1 (552)

Se e assegnato un prezzo iniziale

p0

possiamo trovare ricorsivamente

p1, p2, p3, ......, pn, .....

e possiamo studiare il comportamento di pn

29 marzo 2001 415

La situazione puo essere illustrata dai se-

guenti grafici

29 marzo 2001 416

Le prime iterazioni sono illustrate da

29 marzo 2001 417

e si prosegue con

29 marzo 2001 418

Se supponiamo

• Supponiamo lineare la funzione che defi-

nisce la domanda in 550

f(p) = −ap+ b a, b > 0 (553)

• Supponiamo lineare la funzione che defi-

nisce l’offerta in 551

g(p) = cp+ d c, d > 0 (554)

Possiamo descrivere esplicitamente l’andamen-

to del prezzo.

29 marzo 2001 419

Si ha Dn = −apn + b

Sn+1 = cpn + d

Sn+1 = Dn+1

(555)

possiamo ricavare il prezzo di equilibrio dalla

Sn+1 = Dn+1

Avremo

−apn+1 + b = cpn + d (556)

per cui

pn+1 = −c

apn +

b− d

a= Apn +B (557)

29 marzo 2001 420

Un Modello Classicodi crescita economica

in un’economia in espansione

Possiamo considerare il reddito nazionale co-stituito da due componenti i consumi e gliinvestimenti

Indichiamo con

• Yn Il Reddito

• Cn I Consumi

• In Gli Investimenti

relativi all’anno n

Avremo che

Yn = Cn + In (558)

29 marzo 2001 421

Possiamo supporre che

Cn = c+mYn (559)

Cioe i consumi sono proporzionali al reddi-

to a meno di una quota fissa di consumi

irrinunciabili.

La costante m misura la propensione al con-

sumo mentre c misura il consumo minimale

in assenza di reddito.

Supponiamo

c > 0 , 0 < m < 1

29 marzo 2001 422

L’incremento di reddito e legato agli investi-

menti da

∆Yn = Yn+1 − Yn = rIn (560)

Possiamo ricavare che

Yn+1 − Yn = rIn =

= r(Yn − Cn) =

= rYn − r(c+mYn) (561)

e quindi

Yn+1 = [1 + r(1−m)]Yn − rc (562)

29 marzo 2001 423

Equivalentemente si puo ricavare che

In+1 − In = (Yn+1 − Cn+1)− (Yn − Cn) =

= (Yn+1 − Yn)− (Cn+1 − Cn) (563)

e quindi dalle 559 e 561

In+1 − In =

= rIn −m(Yn+1 − Yn) =

= rIn −mrIn (564)

Infine

In+1 = [(1 + r(1−m)]In (565)

29 marzo 2001 424

L’equazione 565 e una equazione alle diffe-

renze finite di primo ordine e la sua soluzione

puo essere espressa come segue.

In = [(1 + r(1−m)]nI0 (566)

La soluzione dell’equazione 562 e data da

Yn = [(1 + r(1−m)]n(Y0 − Y ∗) + Y ∗ (567)

dove Y ∗ = c1−m

29 marzo 2001 425

Un Modellodi Crescita del Reddito Nazionale

Possiamo considerare il reddito nazionale co-

stituito da tre componenti i consumi gli in-

vestimenti e le spese governative

Indichiamo con

• Yn, il reddito nazionale

• Cn, i consumi

• In, gli investimenti

• Gn,le spese governative.

relativi all’anno n

29 marzo 2001 426

Avremo che

Yn = Cn + In +Gn (568)

Seguendo un classico articolo di Samuelson

possono essere fatte le seguenti ipotesi.

• I consumi in un certo periodo sono pro-

porzionali al reddito del precedente perio-

do.

• gli investimenti in un certo periodo so-

no indotti in maniera proporzionale al-

l’aumento dei consumi in quel periodo

rispetto al periodo precedente.

• Le spese governative sono costanti.

29 marzo 2001 427

Possiamo tradurre le ipotesi fatte mediante

le seguenti equazioni

•

Cn = αYn−1

•

In = β(Cn − Cn−1)

•

Gn = 1

29 marzo 2001 428

Possiamo allora ricavare dalla 568 la seguen-

te equazione

Yn = αYn−1 + β(Cn − Cn−1) + 1

= αYn−1 + β(αYn−1 − αYn−2) + 1 (569)

Yn = α(1 + β)Yn−1 − αβYn−2 + 1 (570)

Si puo verificare che tale successione ha un

andamento oscillante.

29 marzo 2001 429

Andamento del Reddito Nazionale in

corrispondenza

dei valori α = 3 e β = 2 con dati iniziali

I0 = 1 ed I1 = .3

29 marzo 2001 430

Andamento del Reddito Nazionale in

corrispondenza

dei valori α = 1 e β = 5 con dati iniziali

I0 = 1 ed I1 = .3

29 marzo 2001 431

Andamento del Reddito Nazionale in

corrispondenza

dei valori α = 5 e β = 1 con dati iniziali

I0 = 1 ed I1 = .3

29 marzo 2001 432

Andamento del Reddito Nazionale in

corrispondenza

dei valori α = 5 e β = 1

con dati iniziali I0 = 1 ed I1 = 3

29 marzo 2001 433

Un ModelloPer la Gestione delle Scorte

Ogni impresa produce beni per due scopi

1. Per destinarli alla vendita

2. Per destinarli al mantenimento di livelli di

scorta ottimali.

29 marzo 2001 434

Supponendo di usare appropriate unita di mi-

sura possiamo individuare nella produzione

totale tre distinte componenti la parte de-

stinata alla vendita, la parte destinata alle

scorte ed una parte costante non indotta.

Indichiamo con

• Yn, la produzione totale

• Un la parte destinata alla vendita,

• Sn, la parte destinata alle scorte

• V0, la quota di produzione fissa

relative all’anno n

29 marzo 2001 435

Avremo che

Yn = Un + Sn + V0 (571)

Un’estesa analisi dei cicli di mantenimento

delle scorte e stata proposta da L.A.Metzeler.

Per illustrarne i rudimenti possiamo comin-

ciare con il considerare alcune ipotesi

• La quota destinata alla vendita e propor-

zionale all’intera produzione del periodo

precedente.

• La quota destinata alle scorte e decisa

all’inizio di ogni periodo ed e destinata a

mantenere costante il livello delle scorte

29 marzo 2001 436

Possiamo tradurre le ipotesi fatte mediantele seguenti equazioni

•

Un = βYn−1

•

Sn = Un−1 − Un−2 = β(Yn−1 − Yn−2)

Possiamo allora ricavare dalla 571 la seguen-te equazione

Yn = βYn−1 + β(Yn−1 − Yn−2) + V0

= 2βYn−1 − βYn−2 + V0 (572)

Yn = 2βYn−1 − βYn−2 + V0 (573)

Si puo verificare che tale successione ha unandamento oscillante.29 marzo 2001 437

Andamento delle Scorte in corrispondenza

dei valori α = .6 e V0 = 1

con dati iniziali Y0 = 1 ed Y1 = 3

29 marzo 2001 438

Andamento delle Scorte in corrispondenza

dei valori α = .6 e V0 = 1

con dati iniziali Y0 = 3 ed Y1 = 1

29 marzo 2001 439

Andamento delle Scorte in corrispondenza

dei valori α = .6 e V0 = 1

con dati iniziali Y0 = 4 ed Y1 = 3

29 marzo 2001 440

Le equazioni alle differenze

29 marzo 2001 441

Le equazioni alle differenzedel primo ordine

Sia a ∈ R ed fn una successione assegnata,

Consideriamo il problema di determinare una

successione yn in modo che

yn = ayn−1 + fn ∀n ∈ N (574)

Il numero a e la successione fn sono i dati

dell’equazione

Parallelamente

Sia a ∈ R ed fn una successione assegnata,

Consideriamo il problema di determinare una

successione yn in modo cheyn+1 = ayn + fn ∀n ∈ Ny0 = α

(575)

Il numero a e la successione fn sono i dati

dell’equazione

29 marzo 2001 442

Vale il seguente teorema di esistenza ed uni-

cita

Teorema 1 Se a ∈ R e se fn e una suc-

cessione assegnata, esiste una ed una sola

successione yn tale cheyn+1 = ayn + fn ∀n ∈ Ny0 = α

(576)

Dimostrazione

L’esistenza della soluzione segue dal principio

di induzione.

Infatti

1. y0 e definito dal dato iniziale

2. se e definito yk per k ≤ n allora la 574

permette di definire yn+1

29 marzo 2001 443

Per quanto riguarda l’unicita se yn e zn sono

soluzioni della 575 si hayn+1 = ayn + fn ∀n ∈ Ny0 = α

(577)

zn+1 = azn + fn ∀n ∈ Nz0 = α

(578)

e posto un = yn − zn si verifica subito sot-

traendo membro a membro che un soddisfa

le condizioniun+1 = aun ∀n ∈ Nu0 = 0

(579)

Si verifica facilmente per induzione che

un = 0 ∀n ∈ N

e quindi che

yn = zn

29 marzo 2001 444

Cominciamo a determinare le soluzioni del-

l’equazione 574 nel caso in cui

fn = 0

avremo che deve essereyn+1 = ayn ∀n ∈ Ny0 = α

(580)

Quindi

y0 = α (581)

y1 = ay0 = αa (582)

y2 = ay1 = αa2 (583)

y3 = ay2 = αa3 (584)

e evidente che

yn = αan (585)

e si puo verificare per induzione che la 585

risolve effettivamente il problema dato.

29 marzo 2001 445

Consideriamo ora l’equazione completa ad esem-

pio nel caso in cui

fn = bn (586)

Consideriamo cioe il problemayn+1 = ayn + bn ∀n ∈ Ny0 = α

(587)

Si ha

y0 = α

y1 = ay0 + 1 = αa+ 1

y2 = ay1 = αa2 + a+ b

y3 = ay2 = αa3 + a2 + ab+ b2

y4 = αa4 + a3 + a2b+ ab2 + b3 (588)

29 marzo 2001 446

ed e evidente che

yn = αan+an−1+ ban−2+ · · ·+abn−2+ bn−1

= αan +an − bn

a− b=

=(α+

1

a− b

)an +

1

b− abn (589)

e si puo verificare per induzione che la 589

risolve effettivamente il problema dato.

29 marzo 2001 447

La 589

yn =(α+

1

a− b

)an +

1

b− abn (590)

suggerisce alcune considerazioni.

1. Nella soluzione si riconosce

• una parte (α+

1

b− a

)an

che risolve l’equazione omogenea e checontiene una costante α che individuail dato iniziale

• ed una parte

1

b− abn

che risolve l’equazione completa e chenon dipende da nessuna costante.

2. La soluzione non e definita nel caso in cuia = b

29 marzo 2001 448

Il caso in cui b = a, non puo essere tratta-

to in questo modo (b − a = 0), in tal caso

l’equazione diventayn+1 = ayn + an ∀n ∈ Ny0 = α

(591)

e possiamo verificare che si ha

y0 = α (592)

y1 = ay0 + 1 = αa+ 1 (593)

y2 = ay1 = αa2 + a+ a = αa2 + 2a (594)

y3 = ay2 = αa3 + 2a2 + a2 = αa3 + 3a2

(595)

y4 = αa4 + 3a3 + a3 = αa4 + 4a3 (596)

29 marzo 2001 449

Si intuisce e si prova per induzione che

yn = αan + nan−1 (597)

infatti

• y0 = αa0 + 0a−1 = α

• se yn = αan + nan−1 allora

yn+1 = a(αan + nan−1) + an =

αan+1 + (n+ 1)an (598)

29 marzo 2001 450

Studiando qualche equazione alle differenze

del primo ordine abbiamo messo in eviden-

za una struttura vettoriale nello spazio delle

soluzioni.

Piu precisamente le soluzioni formano uno

spazio lineare affine di dimensione 1.

Risultati simili si possono ottenere anche per

equazioni alle differenze di ordine qualsiasi.

Procederemo ad illustrare il caso delle equa-

zioni del secondo ordine; i risultati per le

equazioni di ordine arbitrario sono facilmente

intuibili.

29 marzo 2001 451

Le equazioni alle differenzedel secondo ordine

29 marzo 2001 452

Consideriamo un’equazione alle differenze di

secondo ordine:

Se a, b ∈ R e se fn e una successione assegna-

ta, determinare una successione yn in modo

che

yn+2 = ayn+1 + byn + fn ∀n ∈ N (599)

I numeri a, b e la successione fn sono i dati

dell’equazione

Parallelamente consideriamo il problema di

Se a, b ∈ R e se fn e una successione assegna-

ta, determinare una successione yn in modo

cheyn+2 = ayn+1 + byn + fn ∀n ∈ Ny0 = α

y1 = β

(600)

I numeri a, b, y0, y1 e la successione fn sono i

dati dell’equazione

29 marzo 2001 453

Osserviamo che corre una stretta analogia

tra l’equazione alle differenze finite e le equa-

zioni lineari.

Infatti, se definiamo

∆yn = yn+1 − yn (601)

e

∆2yn = ∆(∆yn) = ∆yn+1 −∆yn = (602)

= yn+2 − 2yn+1 + yn (603)

la 599 puo essere riscritta tenendo conto che

yn+2 − ayn+1 − byn =

= yn+2 − 2yn+1 + yn + (2− a)yn+1 −− (2− a)yn + (2− a− b)yn =

∆2yn + (2− a)∆yn + (3− a− b)yn (604)

e diventa del tipo

∆2yn + a1∆yn + a0yn = fn (605)

29 marzo 2001 454

Possiamo provare che

Teorema 2 Se a, b ∈ R e se fn e una

successione assegnata, Esiste una ed

una sola successione yn tale cheyn+2 = ayn+1 + byn + fn ∀n ∈ Ny0 = α

y1 = β

(606)

29 marzo 2001 455

L’esistenza della soluzione si puo far discen-

dere dal principio di induzione.

Infatti

1. y0, y1 sono definiti

2. se sono definiti yn−1 ed yn−2 allora la 599

permette di definire yn

Per quanto riguarda l’unicita se yn e zn sono

soluzioni della 600 si ha

29 marzo 2001 456

yn+2 = ayn+1 + byn + fn ∀n ∈ Ny0 = α

y1 = β

(607)

zn+2 = azn+1 + bzn + fn ∀n ∈ Nz0 = α

z1 = β

(608)

e possiamo verificare per sottrazione membro

a membro che

un = yn − zn

soddisfa le condizioniun+2 = aun+1 + bun ∀n ∈ Nu0 = 0

u1 = 0

(609)

29 marzo 2001 457

Si verifica facilmente per induzione che

un = 0 ∀n ∈ N

e quindi che

yn = zn

Per caratterizzare le soluzioni di 599 introdu-

ciamo l’equazione omogenea associata, che

sara data da

yn+2 = ayn+1 + byn ∀n ∈ N (610)

29 marzo 2001 458

Definizione 1

S = {yn : yn+2 = ayn+1 + byn} (611)

T = {yn : yn+2 = ayn+1 + byn + fn}(612)

Lo spazio S e uno spazio vettoriale di

dimensione 2

infatti se yn, zn ∈ S avremo che

yn+2 = ayn+1 + byn

zn+2 = azn+1 + bzn (613)

e quindi

αyn + βzn ∈ S

29 marzo 2001 459

Inoltre si vede che se yn, zn ∈ T avremo

yn+2 = ayn+1 + byn + fn

zn+2 = azn+1 + bzn + fn (614)

per cui

(yn+2 − zn+2) = a(yn+1 − zn+1) + b(yn − zn)

(615)

e

yn − zn ∈ S

Ne deduciamo che S e uno spazio vettoriale

e T e uno spazio lineare affine in quanto per

ogni fissata z∗n ∈ T e per ogni yn ∈ T si ha

yn − z∗n ∈ S

e quindi

yn ∈ z∗n + S

29 marzo 2001 460

Sia

Definizione 2 Definiamo

T : S → R2 (616)

mediante la

yn ∈ S 7→ T (yn) = (y0, y1) ∈ R2 (617)

e

Teorema 3 T : S → R2

• e lineare

• e surgettiva

• e iniettiva

Pertanto e un isomorfismo tra spazi vettoriali

e risulta dim S = 2

29 marzo 2001 461

Infatti

• T e lineare in quanto se yn, zn ∈ S, α, β ∈R avremo

T (αyn+βzn) = (αy0+βz0, αy1+βz1) =

= α(y0, y1) + β(z0, z1) =

= αT (yn) + βT (zn) (618)

• T e surgettiva, in quanto se (α, β) ∈ R2,

per il teorema 2 ∃yn tale cheyn+2 = ayn+1 + byn ∀n ∈ Ny0 = α

y1 = β

(619)

Per cui yn ∈ S ed inoltre

T (yn) = (y0, y1) = (α, β)

29 marzo 2001 462

• T e iniettiva in quanto se yn, zn ∈ S e

T (yn) = T (zn) avremo cheyn+2 = ayn+1 + byn ∀n ∈ Ny0 = α

y1 = β

(620)

e zn+2 = azn+1 + bzn ∀n ∈ Nz0 = α

z1 = β

(621)

e per il teorema 2 la soluzione e unica

Pertanto il problema di descrivere lo spazio

S si riduce a trovare due soluzioni y1n ed y2nlinearmente indipendenti dell’equazione

yn+2 = ayn+1 + byn ∀n ∈ N (622)

ed una soluzione della

yn+2 = ayn+1 + byn + fn ∀n ∈ N (623)

29 marzo 2001 463

A questo scopo si cercano soluzioni della

yn+2 = ayn+1 + byn ∀n ∈ N (624)

della forma

yn = ζn = (eα+ i β)n =

= en(α+ i β) = eαn(cosβn+ i sinβn)

ζ ∈ C, α, β ∈ R (625)

Dovra risultare

ζn+2 = aζn+1 + bζn (626)

e, scartando la soluzione banale nulla

ζ2 = aζ + b (627)

29 marzo 2001 464

Si possono presentare tre casi:

1. l’equazione 626 ha due radici reali e di-

stinte λ, µ ∈ R , λ 6= µ;

2. l’equazione 626 ha due radici complesse

e coniugate α± β = ρ(cos θ ± i sin θ)

3. l’equazione 626 ha due radici reali e coin-

cidenti λ ∈ R con molteplicita 2

cui corrispondono altrettante coppie di suc-

cessioni che risolvono l’equazione

29 marzo 2001 465

1. le successioni

yn = λn zn = µn

sono soluzioni dell’equazione

2. Le successioni

yn = ρn cos θn zn = ρnsinθn

sono soluzioni del problema

3. Le successioni

yn = λn zn = nλn

sono soluzioni del problema

29 marzo 2001 466

Nel caso delle radici reali e distinte e evidente

come si ricavano le due soluzioni

yn = λn zn = µn

mentre e meno chiaro come si procede negli

altri due casi.

Consideriamoli piu attentamente.

Sia

α± i β = ρ(cos θ ± i sin θ)

la coppia di soluzioni complesse e coniugate

dell’equazione 626

Avremo allora che le successioni complesse

ξn = (α+ i β)n ζn = (α− i β)n

risolvono l’equazione 574.

29 marzo 2001 467

Infatti

ξn = eαn(cos θn+ sin θn)

ζn = eαn(cos θn− sin θn) (628)

e

ξn = ρn(cos θn+ sin θn)

ζn = ρn(cos θn− sin θn) (629)

Per la linearita anche

yn = ρn cos θn zn = ρnsinθn

soddisfano l’equazione 599

Nel caso in cui l’equazione 626 ammetta una

sola soluzione reale con molteplicita 2 pos-

siamo riscrivere la 599

29 marzo 2001 468

yn+2 = 2ayn+1 − a2yn = 0 (630)

da cui

yn+2 − 2ayn+1 + a2yn = 0 (631)

e

yn+2 − ayn+1 − ayn+1 + a2yn = 0 (632)

(yn+2 − ayn+1)− a(yn+1 − ayn) = 0 (633)

e posto un = yn+1 − ayn

un+1 − aun = 0 (634)

Se ne ricava che

un = kan

yn+1 − ayn = kan

ed ancora

yn = c1an + c2na

n

29 marzo 2001 469

Possiamo affermare che lo spazio T dell’e-

quazione completa puo essere descritto se

• si e in grado di descrivere lo spazio Sdell’equazione omogenea associata

• e noto almeno un elemento y∗ ∈ T , cioe

una soluzione particolare dell’equazione

completa

E quindi importante essere in grado di tro-

vare una soluzione particolare dell’equazione

completa.

29 marzo 2001 470

A questo scopo disponiamo di due strumenti.

Se

fn = ρn (p1(n) sin θn+ p2(n) cos θn)

con p1, p2 polinomi di grado al piu m, pos-

siamo cercare una soluzione particolare nella

forma

y∗n = nµρn (q1(n) sin θn+ q2(n) cos θn)

dove q1, q2 sono polinomi di grado m con

coefficienti da determinarsi e µ e la molte-

plicita di ρ± ı θ nel caso in cui si verifichi che

α± ıβ = ρ± ı θ

29 marzo 2001 471

In caso contrario possiamo procedere in ma-

niera del tutto simile al metodo di variazio-

ne delle costanti arbitrarie di Lagrange come

segue

Se yn ed zn sono soluzioni dell’equazione da-

ta, possiamo cercare una soluzione dell’equa-

zione completa nella forma

y∗n = λnyn + µnzn

imponendo opportune condizioni su λk e µk.

Poiche yn e zn solo soluzioni dell’equazione

omogenea

∆2yn + a1∆yn + a0yn = 0

∆2zn + a1∆zn + a0zn = 0

29 marzo 2001 472

Possiamo anche calcolare che

∆y∗n = λn∆yn + µn∆zn + ∆λnyn+1 + ∆µnzn+1

(635)

e, posto

∆λnyn+1 + ∆µnzn+1 = 0

∆2y∗n = λn∆2yn + µn∆

2zn +

+ ∆λn∆yn+1 + ∆µn∆zn+1 (636)

Onde, per sostituzione nell’equazione data,

tenendo conto che yn e zn sono soluzioni

dell’omogenea, si ha∆yn+1∆λn + ∆zn+1∆µn = bn

yn+1∆λn + zn+1∆µn = 0(637)

ed e possibile da esse ricavare

∆λn e ∆µn

e λn e µn successivamente.

29 marzo 2001 473

E anche importante studiare la lineare indi-

pendenza di due (o piu successioni).

A questo scopo ricordiamo che yn e zn sono

soluzioni linearmente indipendenti dell’equa-

zione che stiamo considerando se

λyn + µzn = 0 ∀n ∈ N ⇐⇒ λ = µ = 0

(638)

Cio implica che il sistemaλyn + µzn = 0

λyn+1 + µzn+1 = 0(639)

oppure il sistema

λyn + µzn = 0

λ(yn+1 − yn) + µ(zn+1 − zn) = 0(640)

ha la sola soluzione nulla

29 marzo 2001 474

Ovviamente il sistema 640 puo essere riscrit-

to nella formaλyn + µzn = 0

λ∆yn) + µ∆zn = 0(641)

ed ammette la sola soluzione nulla.

Ne viene che

det =

(yn zn = 0yn+1 zn+1

)6= 0(

det

(yn zn

∆yn ∆zn

)6= 0

)(642)

per ogni n ∈ N

29 marzo 2001 475

Pertanto si ha anche che

det =

(y0 z0y1 z1

)6= 0

(det

(y0 z0

∆y0 ∆z0

)6= 0

)(643)

D’altro canto se quest’ultima condizione e

verificata. allora (y0, y1), (z0, z1) ∈ R2 sono

linearmente indipendenti e tali che

T (yn) = (y0, y1) e T (zn) = (z0, z1)

e quindi

yn = T−1(y0, y1) e zn = T−1(z0, z1)

e ed essendo T un isomorfismo e (y0, y1), (z0, z1)

linearmente indipendenti anche yn e zn sono

linearmente indipendenti.

La matrice(y0 z0y1 z1

)oppure

(y0 z0

∆y0 ∆z0

)(644)

prende il nome di matrice di Casorati e pos-

siamo riassumere quanto detto come segue.

29 marzo 2001 476

Siano yn, zn due soluzioni della 599.

Sono condizioni equivalenti

• yn, zn sono linearmente indipendenti

•

det =

(yn zn = 0yn+1 zn+1

)6= 0

det

((yn zn

∆yn ∆zn

))6= 0 (645)

per ogni n ∈ N•

det =

(y0 z0y1 z1

)6= 0

det

((y0 z0

∆y0 ∆z0

))6= 0 (646)

29 marzo 2001 477

Un esempio interessante e costituito dalla

successioneyn+2 = 2yn+1 cos 1− yn

y0 = 0

y1 = sin1

(647)

Si ha

P (λ) = λ2 − 2λ cos 1 + 1 = 0

e

λ = cos1 + i sin 1 = ei

da cui, tenendo anche conto delle condizioni

iniziali

yn = sinn.

29 marzo 2001 478

Qualche Profilo Storico

29 marzo 2001 479

Francois Edouard Anatole Lucasnato ad Amiens il 4 Aprile 1842

morto a Parigi il 3 Ottobre 1891

29 marzo 2001 480

Lucas fu educato alla scuola Normale Supe-

riore di Amiens e lavoro all’Osservatorio di

Parigi

Fu Ufficiale di artiglieria durante la guerra

Franco-Prussiana (1870-1871) e successiva-

mente insegno nei licei a Parigi.

si occupo principalmente di Teoria dei Nume-

ri e dimostro che 2127−1 e un numero primo

usando metodi da lui sviluppati.

2127 − 1 e a tutt’oggi il piu grande numero

primo trovato senza l’ausilio di un computer.

Lucas e anche ben noto per l’invenzione della

Torre di Hanoi e di altri giochi di carattere

matematico.

29 marzo 2001 481

Leonardo da Pisa detto Fibonacci

cioe

figlio di Bonaccio

29 marzo 2001 482

• Nacque a Pisa attorno al 1170

• Morı a Pisa attorno al 1250

• Fu educato in Nord Africa da precettori

mussulmani

• Ebbe modo di conoscere ed apprezzare il

sistema di numerazione indo-arabica che

introdusse per primo in Europa.

• Le sue opere maggiori sono

– Liber Abaci (1202)

– Practica Geometriae (1220)

– Liber Quadratorum (1225)

29 marzo 2001 483

Paul Anthony Samuelson

Nacque a Gary nell.Indiana nel 1915

Premio Nobel per l’economia nel 1970

29 marzo 2001 484

Il Suo lavoro piu noto e contenuto nel libro

Foundation of Economic Analysis

Ha incrementato in maniera notevole l’uso

della Matematica nell’economia

Secondo Samuelson quasi ogni comportamen-

to puo essere interpretato nei termini di un

problema di massimo o di minimo soggetto

a vincoli.

29 marzo 2001 485

Thomas Malthus (1766-1834) era un eco-nomista inglese che diede inizio allo studiodemografico organizzato. Il suo contributoprincipale e stato appunto il modello di cre-scita delle popolazioni che pubblico nel suolibro ”An Essay on the Principle of Polula-tion” (1798). Nella sua opera egli sostenneche la poverta e causata dal fatto che la po-polazione cresce piu velocemente dei mezzidi sussistenza. Come controlli sulla crescitadella popolazione Malthus accetto soltanto laguerra, la carestia e le malattie, ma in seguitoaggiunse anche la condotta morale

29 marzo 2001 486

Pierre Francois Verhulst, era un matematico

Belga vissuto tra il 1804 ed il 1849

Lavoro sulla teoria dei numeri e sulla statisti-

ca, in particolare sulla statistica sociale.

Le sue ricerche sulla crescita delle popolazio-

ni sono il suo risultato piu rappresentativo.

Basandosi sulla sua teoria predisse che nel

1994 la popolazione del Belgio sarebbe stata

di 9.400.000 individui.

In realta ne risultarono 10.228.000, ma pos-

siamo dire che le le sue previsioni sono di

tutto rispetto se teniamo conto del fenomeno

dell’immigrazione.

29 marzo 2001 487