Problemi di Controllo ad orizzonte infinito

Problemi di Controllo ad orizzonte
infinito-soluzioni numeriche e applicazione ad un
modello di crescita economica sostenibile
Docente:Alessandra Cutrı̀
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Ricostruzione mediante procedimenti numerici della
funzione valore
Ricostruire la funzione valore di un problema di controllo ad
orizzonte infinito
Versione discreta del principio di programmazione dinamica
Equazione di Bellman discreta
Metodo delle approssimazioni successive e Metodo alle
differenze finite
basati su Teorema delle Contrazioni
algoritmo convergente alla funzione valore V
Testi consigliati:
M.Bardi-I.Capuzzo Dolcetta: “Optimal control and viscosity
solutions of Hamilton-Jacobi-Bellman equations.”, Systems &
Control: Foundations & Applications. Birkhäuser Boston Inc.,
Boston, MA, 1997.
I.Capuzzo Dolcetta-M.Falcone: “Discrete dynamic
programming and viscosity solutions of the Bellman
equation”, Annales de l’I.H.P. (1989).
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Riprendiamo il Problema di controllo ad orizzonte infinito
Problema ad orizzonte infinito: minimizzare il funzionale costo
definito su RN × A
Z ∞
J(x, α) =
l(y (t; x, α), α(t))e −λt dt
0
dove l : RN × A → [0, +∞) è il costo corrente che verifica:
l continua nei suoi argomenti
∀a ∈ A, ∀x, z ∈ RN
|l(x, a) − l(z, a)| ≤ Ll |x − z|
klk∞ ≤ Ml
e λ > 0 è il tasso di sconto (che corrisponde al fatto che costi
futuri pesano meno).
la funzione valore è definita da
V (x) = inf J(x, α)
α∈A
A. Cutrı̀
11-01-2016, 13-01-2016
x ∈ RN
Metodi Matematici per l’ingegneria–Ing. Gestionale
Ricordiamo Sistema dinamico controllato
Le traiettorie sono le soluzioni del sistema dinamico controllato in
RN :
ẏ (t; x, α) = f (y (t; x, α), α(t)) t > 0
(1)
y (0) = x
(a volte indicheremo con y (t, α) o anche y (t) la soluzione) dove
f : RN × A → RN , A è un compatto di RN ed il controllo
α : [0, ∞] → A è una funzione misurabile
la dinamica f soddisfa le ipotesi che garantiscono esistenza ed
unicità della soluzione globale di (1) (per esempio f è limitata e
Lipschitziana in x), cioè:
f continua, ∃Mf > 0 , tale che , |f (x, a)| ≤ Mf , ∀x ∈ RN , ∀a ∈ A
∃Lf > 0 tale che |f (x, a)−f (z, a)| ≤ Lf |x−z|
A. Cutrı̀
11-01-2016, 13-01-2016
∀x, z ∈ RN , ∀a ∈ A
Metodi Matematici per l’ingegneria–Ing. Gestionale
Principio di Programmazione dinamica e soluzione di
viscosità dell’equazione H-J
Come determinare V e soprattutto la strategia ottima?Attraverso il
principio di programmazione dinamica:
Principio di Programmazione dinamica (PPD): Sotto le ipotesi
indicate su f e su l, per ogni fissato T > 0, V soddisfa:
Z
V (x) = inf {
α∈A
0
T
l(y (t; x, α), α(t))e −λt dt+V (y (T , x, α))e −λT }
e si può provare che V è l’unica soluzione viscosità dell’equazione
di Hamilton-Jacobi:
λV (x) + sup{−f (x, a) · DV (x) − l(x, a)} = 0
a∈A
A. Cutrı̀
11-01-2016, 13-01-2016
x ∈ RN
Metodi Matematici per l’ingegneria–Ing. Gestionale
Versione discreta del problema di controllo a orizzonte
infinito
Vediamo come il metodo della programmazione dinamica discreta
sia utile per risolvere il problema continuo.
Discretizzazione temporale: Sia h > 0 un passo temporale
fissato e supponiamo che l’evoluzione del sistema dinamico sia
osservata solo ad una successione di istanti tj = jh,
j = 0, 1, 2, . . . . Supponiamo che nell’intervallo [tj , tj+1 [ la
dinamica ed il costo corrente rimangano costanti, cioè:
f (y (t), α(t)) = f (yj , aj )
t ∈ [tj , tj+1 [
l(y (t), α(t)) = l(yj , aj )
t ∈ [tj , tj+1 [
dove aj = α(tj ) e yj = y (tj ) è definita per ricorrenza da:
y0 = x
yj+1 = yj + hf (yj , aj )
j = 0, 1, 2, . . .
A. Cutrı̀
11-01-2016, 13-01-2016
(2)
Metodi Matematici per l’ingegneria–Ing. Gestionale
La successione yj definita in (2) rappresenta una versione discreta
delle traiettorie del sistema dinamico controllato (1):
ẏ (t; x, α) = f (y (t; x, α), α(t)) t > 0
y (0) = x
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Funzionale costo e funzione valore discreti
funzionale costo associato all’evoluzione discreta (2)
Jh (x, α) = h
∞
X
β j l(yj , aj )
j=0
con β = 1 − λh
funzione valore approssimata
Vh = inf {Jh (x, α)}
α∈Ah
dove Ah ⊂ A è costituito dai controlli costanti in [tj , tj+1 [.
Principio di programmazione dinamica discreto (DPPh ):
∀x = y0 ∈ RN , p = 0, 1, 2, . . .
Vh (x) = inf {h
α∈Ah
A. Cutrı̀
p−1
X
β j l(yj , aj ) + β p Vh (yp )}
j=0
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Equazione soddisfatta da Vh in conseguenza di (DPPh )
Per il (DPPh ), la funzione valore Vh è l’unica soluzione limitata
della seguente equazione funzionale di Bellman
Vh (x) + sup[−βVh (x + hf (x, a)) − hl(x, a)] = 0 x ∈ RN (3)
a∈A
per h ∈ (0, λ1 ). (Basta scegliere p = 1 in (DPPh ))
Se λ > Lf , per h → 0 Vh converge localmente uniformemente
alla funzione valore V del problema continuo ad orizzonte
infinito ed esiste C > 0 tale che
√
kVh − V k∞ ≤ C h
(sotto opportune ipotesi sulla dinamica f e costo corrente l
l’errore diventa di ordine h).
per determinare Vh occorre ridurre il problema (3) ad un
problema finito dimensionale.
⇒ necessaria una discretizzazione della variabile di stato x e
Teorema delle contrazioni (Teorema di punto fisso).
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Sintesi dei controlli ottimi per Vh e V
Supponiamo Vh semicontinua inferiormente allora, fissato x ∈ RN
il sup in (3):
Vh (x) + sup[−βVh (x + hf (x, a)) − hl(x, a)] = 0
a∈A
è raggiunto per ah∗ = ah∗ (x) ∈ A. Definiamo
x0∗ = x
∗
xj+1
= xj∗ + hf (xj∗ , ah∗ (xj∗ )) j = 0, 1, . . .
e poniamo
ah∗ (t) = ah∗ (xj∗ )
t ∈ [tj , tj+1 [
(4)
Allora Se h ∈ (0, λ1 ) e A compatto , il controllo costante a tratti
definito in (4) è ottimo per (3):
Vh (x) = Jh (x, ah∗ )
∀x ∈ RN
e per h → 0, ah∗ converge al controllo ottimo per il problema ad
orizzonte infinito con costo J
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Discretizzazione nella variabile di stato
Per ridurre (3) ad un problema finito dimensionale, è necessaria
una discretizzazione della variabile di stato x. A tale proposito,
supponiamo esista un sottoinsieme aperto e limitato Ω ⊂ RN che
sia invariante per la dinamica f . Costruiamo una triangolazione di
Ω in un numero finito P di simplessi Sj in modo che l’insieme
Ωk = ∪Pj=1 Sj , dove k = maxj diam(Sj ), sia invariante rispetto alle
traiettorie discretizzate:
∃h > 0 : x + hf (x, a) ∈ Ωk
∀(x, a) ∈ Ωk × A
Indichiamo con xi per i = 1, . . . , N i nodi della triangolazione di Ω
e sostituiamo (3) con il sistema di N equazioni
u(xi ) + sup[−βu(xi + hf (xi , a)) − hl(xi , a)] = 0,
i = 1, . . . , N
a∈A
(5)
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Teorema di Banach-Caccioppoli delle Contrazioni
La soluzione di (5) verrà determinata come il punto fisso di
un’opportuna mappa. L’esistenza di tale punto fisso è garantita dal
seguente teorema, molto utile in Analisi Matematica e nelle
applicazioni
Theorem
Sia (X , d) uno spazio metrico completo e T : X → X una
contrazione, cioè esista M ∈ (0, 1) tale che
d(T (x), T (y )) ≤ Md(x, y )
Allora esiste un unico x ∈ X tale che T (x) = x (punto fisso)
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Dimostrazione teorema delle contrazioni
Scegliamo arbitrariamente x0 ∈ X e definiamo la successione per
ricorrenza:
xk+1 = T (xk )
{xk } è una successione di Cauchy in (X , d).
Infatti:
d(x1 , x2 ) = d(T (x0 ), T (x1 )) ≤ Md(x0 , x1 ) = Md(x0 , T (x0 ))
d(x2 , x3 ) = d(T (x1 ), T (x2 )) ≤ Md(x1 , x2 ) = M 2 d(x0 , T (x0 ))
d(xk , xk+1 ) = d(T (xk−1 ), T (xk )) ≤ Md(xk−1 , xk ) = M K d(x0 , T (x0 ))
Per la disuguaglianza triangolare
d(xk , xk+p ) ≤ d(xk , xk+1 ) + · · · + d(xk+p−1 , xk+p )
k −M k+p
d(x0 , T (x0 ))
≤ (M k + M k+1 + · · · + M k+p−1 )d(x0 , T (x0 )) = M 1−M
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Essendo 0 < M < 1 abbiamo
d(xk , xk+p ) ≤
Mk
d(x0 , T (x0 ))
1−M
Poiché M k → 0 per k → ∞, {xk } è di Cauchy
Poiché (X , d) è completo, esiste x ∈ X tale che xk → x in X .
Passando al limite nella relazione xk+1 = T (xk ) ed usando la
continuità di T otteniamo
x = T (x)
il punto fisso è unico: Per assurdo supponiamo esista un altro
punto y ∈ X tale che y = T (y ), allora
d(x, y ) = d(T (x), T (y )) ≤ Md(x, y )
assurdo se d(x, y ) 6= 0 in quanto M < 1
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Soluzione dell’equazione di Bellman Discretizzata nella
variabile di stato
Per ridurre (3) ad un problema finito dimensionale, abbiamo
operato una discretizzazione della variabile di stato x, sotto le
ipotesi di esistenza di un dominio Ω invariante per f . Con una
triangolazione di Ω in un numero finito P di simplessi Sj in modo
che l’insieme Ωk = ∪Pj=1 Sj , dove k = maxj diam(Sj ), sia invariante
rispetto alle traiettorie discretizzate:
∃h > 0 : x + hf (x, a) ∈ Ωk
∀(x, a) ∈ Ωk × A
abbiamo ottenuto da (3) il sistema di N equazioni (5):
u(xi ) + sup[−βu(xi + hf (xi , a)) − hl(xi , a)] = 0,
i = 1, . . . , N
a∈A
dove xi per i = 1, . . . , N sono i nodi della triangolazione di Ω
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Metodo classico delle approssimazioni successive
approssimazione iniziale u(xi ) = ui0 , i = 1, . . . , N
si determina u 0 (xi + hf (xi , a)) mediante interpolazione sui
nodi della griglia
Si definisce la successione ricorsiva:
u n+1 (xi ) = Th (u n ) = inf {βu n (xi +hf (xi , a))+hl(xi , a)}, i = 1, . . . , N
a∈A
Th è una mappa contrattiva, cioè soddisfa
kTh (u) − Th (v )k∞ ≤ βku − v k∞
essendo β < 1
lemma delle contrazioni ⇒ esiste un unico punto fisso per Th
cioè un’unica funzione u che soddisfa
u = Th (u)
⇒ u n converge per n → ∞ ad una soluzione di (5)
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
questo metodo, pur convergendo a partire da qualsiasi scelta di u 0 ,
ha molti difetti:
la convergenza è lenta poiché per h → 0, il coefficiente della
contrazione β → 1
costo computazione elevatissimo: poiché per ottenere
u n+1 (xi ) (quindi ad ogni iterazione) se si hanno a disposizione
M controlli, poiché xi + hf (xi , a) ∈ Ωk è necessario prendere
gli N valori di u n nei nodi e confrontarli con gli M differenti
controlli: perciò necessari N × M confronti ad ogni iterazione
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Metodo alle differenze finite
Si cerca una soluzione di (5)
u(xi ) + sup[−βu(xi + hf (xi , a)) − hl(xi , a)] = 0,
i = 1, . . . , N
a∈A
nello spazio delle funzioni lineari a tratti in Ωk , (xi sono i nodi
della triangolazione)
W = {w ∈ C (Ωk ) : Dw (x) = cj , x ∈ Sj }
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Vale il seguente teorema
Theorem
Sotto le ipotesi standard su f e l, per ogni h ∈ (0, λ1 ) che soddisfa
x + hf (x, a) ∈ Ωk quando (x, a) ∈ Ωk × A, esiste una ed una sola
soluzione vhk ∈ W di (5).
Idea della dimostrazione:
Per ogni funzione affine u k , (5) si scrive
u k (xi ) = inf {β
a∈A
L
X
λj (xi , a)u k (xj ) + hl(xi , a)}, i = 1, . . . , N
j=1
Infatti, essendo Ωk unione di simplessi,
xi + hf (xi , a) =
L
X
λj (xi , a)xj
j=1
0 ≤ λj ≤ 1,
L
X
λj (xi , a) = 1
j=1
uk
Ogni
∈ W , si scrive
PLcome
k
u (xi + hf (xi , a)) = j=1 λj (xi , a)u k (xj )
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
per u affine, (5) è equivalente al problema di punto fisso a
dimensione finita:
U = Th (U)
dove T : RN → RN è definita da
(Th (U))i = min{βΛ(a)U + hL(a)}i , i = 1, . . . , N
a∈A
dove Ui = u(xi ), Λ(a) è una matrice N × N positiva
Λij (a) = λj (xi , a), Li (a) = l(xi , a) per i = 1, . . . , N
Th è una contrazione in RN poiché soddisfa
kTh (U) − Th (V )k ≤ βkU − V k
lemma delle contrazioni ⇒
∃!V ∗
∈
∀U, V ∈ RN ,
RN
β ∈ (0, 1)
tale che Th (V ∗ ) = V ∗
V ∗ = (vhk (x1 ), vhk (x2 ), . . . , vhk (xN ))
interpolando V ∗ otteniamo una soluzione vhk ∈ W di (5)
Si prova poi chevhk converge localmente uniformemente alla
soluzione Vh di (3) quando il diametro k della triangolazione tende
a zero
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Velocità di convergenza
Vale la stima di convergenza:
kvhk (x) − Vh (x)k∞ ≤
Ll
k
λ(λ − Lf ) h
Tenendo conto che Vh converge localmente uniformemente a √
V
soluzione viscosità dell’equazione di H-J continua, con ordine h,
nell’insieme Ω invariante per la dinamica f , se λ > Lf ed f , l
soddisfano le ipotesi standard, si ha per h ∈ (0, λ1 ):
√
∃C > 0 : max |vhk (x) − V (x)| ≤ C h +
x∈Ω
Ll
k
λ(λ − Lf ) h
dove h è il passo di discretizzazione temporale (delle traiettorie) e
k = maxj (diamSj ) è il passo di discretizzazione spaziale
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
U n = Th (U n−1 ) converge qualunque sia la scelta di U 0
convergenza lenta poiché β → 1 quando h → 0,
si utilizza allora procedura di accelerazione monotona: Si
sceglie opportunamente l’approssimazione iniziale U 0
nell’insieme
U ≡ {U ∈ RN : U ≤ Th (U)}
U ≤ V ⇒ Th (U) ≤ Th (V ) implica U n = Th (U n−1 ) è una
successione Monotona crescente e U n ∈ U se U 0 ∈ U
U è chiuso e convesso ed il punto fisso V ∗ di Th rappresenta
l’elemento massimale di U cioè il limite di tutte le successioni
monotone crescenti (anche non strettamente) che partono da
U
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Algoritmo accelerato:
PASSO 0:Si prende U 0 ∈ RN
n := 0
1
PASSO 1: Si calcola U n+ 2 = Th (U n )
1
PASSO 2: Si calcola U n+1 = U n + π(U n+ 2 − U n ) dove
1
π = max{γ ∈ R+ : U n + γ(U n+ 2 − U n ) ∈ U}
PASSO 3: Criterio di arresto: kU n+1 − U n k∞ < ε?
Se SI’ ⇒ STOP
Se NO ⇒ si sostiuisce U n con U n+1 e si torna al PASSO 1
In questo caso la mappa Th viene usata solo per trovare una
direzione di spostamento e la velocità di convergenza monotona a
V ∗ non è strettamente legata a quanto β è vicina a 1
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
tion on the optimal trajectories one can compute the value function. The viscosity
solution in ~t is
Esempio v(x)=
numerico
f (x+l) forx<0,
[3(1-x)
elsewhere.
The solution has one sharp kink at the point x = 0. It can be seen t h a t the numerical solution is close to the exact solution also at that point for a discretization
based on 201 grid points (Figure 2a). Note t h a t at x = 0 one has two different
optimal choices for the control a* = 1 and a* = - 1 and this is the main cause for
the j u m p in the derivative of the solution. Using a higher order method for the
discretization of the system of differential equations and of the cost functional,
the approximate solution can even be improved. In fact, coupling a R u n g e - K u t t a
scheme with a trapezoid rule produces an accurate approximation with just 51 grid
points (see Figure 2b). Several other couplings between a one-step method for ordinary differential equations and a quadrature formula for the cost are possible
(see the hints at the end of this Appendix).
<~
Consideriamo il seguente esempio numerico in R2 (cfr. libro
M.Bardi-I.Capuzzo Dolcetta pag. 483):
Prendiamo il problema di controllo ad orizzonte infinito con
A = B(0, 1) ⊂ R2 , λ = 1, f (x1 , x2 , a) = (a2 x2 , 0),
Ω = (−1, 1) × (−1, 1), l(x1 , x2 , a) = (|x1 | − 1)2 . La funzione valore
TEST 3. Let us consider an example in R 2. We set ~ = ] - 1 , 1[2,
numericaA con
un errore
ε = 0.025 è accurata anche nell’origine
= / ~ ( 0 , 1) c ~2.
A= 1
fl(x,y.a)=ay, una
f2(x,y,a)
=0.
dove presenta
cuspide
FIGURE 3: Value function (h = 0.05, k = 0.025).
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Modello di crescita economica sostenibile
Legame tra processi di crescita economica e qualità dell’ambiente:
È Possibile avere crescita economica e abbattimento dei livelli di
inquinamento?
Vedremo che sviluppo economico e salvaguardia dell’ambiente sono
compatibili!
Protocollo di Kyöto (11-12-1997) entrato in vigore nel 2005
impegnò i paesi industrializzati a ridurre le emissioni di elementi
inquinanti entro il 2012 in misura non inferiore al 5,2% di quelle
registrate nel 1990 prevedendo dei Meccanismi di Flessibilità
(interventi sul mercato adatti a favorire la tutela dell’ambiente
preservando la crescita economica) utili all’acquisizione di crediti
per raggiungimento obiettivo di riduzione emissione
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Meccanismi di Flessibilità nel protocollo di Kyöto
Commercio di emissioni: scambio di crediti di emissione tra
paesi industrializzati, che abbiano ridotto le emissioni di gas
serra più di quanto avrebbero dovuto e Paesi in via di sviluppo
che non siano stati in grado di raggiungere gli obiettivi di
riduzione
Implementazione Congiunta: possibilità per i paesi
industrializzati e per quelli ad economia in transizione di
realizzare progetti per la riduzione di emissioni di gas-serra in
altri paesi dello stesso gruppo e di utilizzare i crediti derivanti
congiuntamente al paese ospite
Meccanismo di sviluppo pulito: come l’implementazione
congiunta ma con realizzazione di progetti in paesi in via di
sviluppo ottenendo in cambio crediti
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Nel 2008 è stato ratificato dall’Unione Europea il pacchetto
20 − 20 − 20: raggiungimento entro il 2020 del 20% della
produzione energetica da fonti rinnovabili, il miglioramento del
20% dell’efficienza ed un taglio del 20% nelle emissioni di
anidride carbonica.
Conferenza sul clima Parigi 2015: contenimento del
riscaldamento globale del pianeta entro i 2◦ C (perseguendo
l’obiettivo ideale di 1, 5◦ C ). Per ottenere il primo risultato si
dovrebbero tagliare le emissioni del 40 − 70% entro il 2050
rispetto a quelle del 2010 (mentre per arrivare a 1, 5◦ le
emissioni andrebbero tagliate tra il 70% ed il 95%)
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Curva ambientale di Kuznets (EKC)
Le curve (EKC) sono curve a campana che descrivono relazioni
EMPIRICHE tra crescita ed inquinamento ambientale. Sono
funzioni che mettono in relazione il PIL pro capite con lo stato
dell’ambiente (per esempio emissione di qualche inquinante).
Secondo queste curve, sperimentali, quando il PIL pro-capite è
basso, scarse risorse vengono destinate alla protezione ambientale e
la crescita avviene a discapito dell’ambiente, ma superato un certo
livello di PIL pro-capite la relazione si inverte e aumento PIL e
protezione ambientale vanno nella stessa direzione
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
La programmazione dinamica applicata ed un modello di crescita
endogena, fornisce una base teorica alla curva EKC ed individua la
percentuale minima di risorse da destinare al progresso tecnologico
per avere una crescita economica sostenibile (con un abbattimento
dell’inquinamento ambientale).
Tesi di Laurea: Dott. Chiara Spada a.a. 2006/2007: ”Analisi dei
modelli di crescita dell’economia Ambientale attraverso la
programmazione dinamica”
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Modellizzazione dell’inquinamento ambientale
Sviluppare modello che unifichi il processo di crescita economica e
l’ambiente
Modulo per descrivere
le caratteristiche del
problema economico
!
+
Modulo ambientale (descrive il
processo di accumulo di
sostanze inquinanti)
L’inquinamento ambientale è un sottoprodotto della
produzione e dei processi di consumo che si sviluppano nel
modulo economico
Le emissioni generate nel modulo economico incidono sul
flusso o sull’accumulo di inquinanti nell’ambiente
L’inquinamento ha effetti deleteri sull’utilità individuale
L’inquinamento può influenzare negativamente la produttività
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Il flusso di emissioni per unità di tempo è legato alla funzione di
produzione Y
Y dipende dal capitale K e dal lavoro effettivo AL
Il capitale K = Ky + Ka dove Ky è il capitale produttivo (che
genera inquinamento) e Ka è il capitale di abbattimento (che
è destinato alla riduzione di sostanze inquinanti attraverso lo
sviluppo di nuove tecnologie e conoscenze)
Y = F (Ka , Ky , AL)
il flusso di emissioni Z = ϕ(Ka )Y dove ϕ è il coefficiente di
emissione unitaria (cioè per unità di prodotto)
Il livello di sostanze inquinanti accumulate nell’ambiente (gas serra,
metalli pesanti etc.) P segue una legge di transizione:
Ṗ = Z − mP +
P2
1 + P2
dove m è il tasso di decadimento esponenziale dell’inquinamento e
P2
è un termine di correzione non lineare.
1+P 2
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Modelli di crescita economica
Esistono due classi di modelli di crescita economica:
Modelli di crescita esogena: il progresso tecnologico è una
variabile esogena (cioè il progresso tecnologico non è spiegato
dal modello di crescita) e la disutilità dovuta all’inquinamento
non è tenuta in considerazione (per esempio nel modello di
Solov). In questo caso non è possibile prevedere crescita
economica sostenibile (senza accumulo di inquinamento)
Modelli di crescita endogena: Il motore della crescita è
l’acquisizione di nuove conoscenze ed il progresso tecnologico
è endogeno. I tassi di crescita possono essere influenzati da
politiche governative come tassazione, regolazione dei mercati
internazionali e possono rimanere positivi se la produttività del
capitale (definito in senso generalizzato in modo da includere
anche il capitale umano) non tende a zero nel lungo periodo o
se la produzione di conoscenza è caratterizzata da guadagni
crescenti. Di questa classe di modelli fa parte il modello AK
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Modello AK con abbattimento e feedback non lineare
La funzione di produzione aggregata (con popolazione costante ed
in assenza di cambiamenti tecnologici esogeni)
Y = AK
dove K tiene conto anche del capitale umano ed A > 0 è un
indice di sviluppo che tiene conto anche del progresso tecnologico
(può essere identificato con la crescita annua di un paese). Il
problema di pianificazione sociale può essere formalizzato come:
Z ∞
max
e −λt U(c, P)dt
c(t)
0
dove c(t) rappresentano i consumi (variabile di controllo), U
rappresenta la funzione di utilità che dipende dai consumi ma in
cui si introduce una disutilità dovuta all’inquinamento.
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
La funzione di utilità
U(c, P) =
cσ
Pγ
−
,
σ
γ
σ ∈ (0, 1) , γ > 1
σ rappresenta la propensione dell’individuo a distribuire i
consumi nel tempo. Più grande è σ, maggiore è la necessita di
consumare un bene nel presente.
γ rappresenta l’ importanza che il consumatore attribuisce al
livello di inquinamento nella valutazione del suo benessere
(maggiore è la disutilià associata all’inquinamento, più γ è
grande
(cfr. e.g. A.Xepapadeas: Economic Growth and the Environment,
2003)
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
La dinamica:
(
K̇ = AKy − c − δ(Ky + Ka )
P2
Ṗ = ϕKy − ψKa − mP + 1+P
2
ϕ: percentuale di capitale produttivo (Ky ) che produce
inquinamento (verosimile scelta ϕ = 0, 1)
ψ: percentuale di capitale di abbattimento (Ka ) che provoca
riduzione delle emissioni inquinanti (verosimile scelta ψ = 0, 5)
m: tasso di decadimento esponenziale dell’inquinamento
(verosimile scelta m = 0, 3)
A = 1.02 (crescita annua di un paese)
δ: tasso di deprezzamento del capitale (verosimile scelta
δ = 0, 2)
Poiché K = Ka + Ky ⇒ Ka = ηK , Ky = (1 − η)K con
η ∈ [0, 1] che rappresenta la percentuale di capitale investito
in Ricerca&Sviluppo. La variabile η è la più importante al fine
di ottenere sviluppo sostenibile.
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
quindi il sistema dinamico controllato diventa:


 K̇ (t) = A(1 − η − δ)K (t) − c(t)
Ṗ(t) = [ϕ(1 − η) − ψη]K (t) − mP(t) +


(K (0), P(0)) = x > 0
P 2 (t)
1+P 2 (t)
(6)
ed il controllo-consumo c(t) ∈ [0, AK (t)]. Il funzionale benessere
Z ∞
P(t)γ
c(t)σ
−
]dt
J(x, c) =
e −λt [
σ
γ
0
e la funzione valore
V (x) =
sup
J(x, c)
{0≤c(t)≤cmax }
risolve l’equazione di H-J:
λV (x) −
A. Cutrı̀
sup
c∈[0,cmax ]
{DV (x) · f (x, c) + [
11-01-2016, 13-01-2016
Pγ
cσ
−
]} = 0
σ
γ
Metodi Matematici per l’ingegneria–Ing. Gestionale
Definiamo W (x) = −V (x), allora W è la funzione valore di un
problema di minimo di un funzionale costo (−J) di tipo orizzonte
infinito e si può applicare il metodo numerico introdotto all’inizio
di queste slides per determinare la soluzione. In particolare, W è
l’unica soluzione di viscosità di
λW (x) +
sup
c∈[0,cmax ]
{−DW (x) · f (x, c) + [
Pγ
cσ
−
]} = 0
σ
γ
Tesi di Laurea: Dott. Chiara Spada a.a. 2006/2007: ”Analisi dei
modelli di crescita dell’economia Ambientale attraverso la
programmazione dinamica”
utilizzata libreria HJPAC (by P.Lanucara, M.Rorro, CASPUR)
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Esperimenti numerici
Vediamo alcuni risultati numerici ottenuti utilizzando la libreria
HJPACK basata sugli algoritmi precedentemente descritti.
Verranno confrontati
funzione valore e controllo ottimo (consumo)
andamento del capitale in funzione del tempo (linea blu)
andamento dell’inquinamento in funzione del tempo (linea
verde)
curva di Kuznets: inquinamento in funzione del PIL pro-capite
facendo variare i parametri:
η (quota capitale destinata a sviluppo di tecnologie pulite)
σ (propensione a distribuire i consumi nel tempo)
γ (disutilità associata a inquinamento)
OSS: σ e γ influenzano la politica dei consumi e la funzione valore
ma non l’andamento del capitale o dell’inquinamento in funzione
del tempo né la curva di Kuznets (come è ragionevole, essendo
parametri soggettivi).
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale
Economia sostenibile se η = 0, 2
Conclusioni: La scelta η = 0, 2 (cioè destinare il 20% del capitale
allo sviluppo di tecnologie pulite) produce
politica di consumi ottimale non troppo restrittiva
crescita incessante del capitale anche a lungo termine
inquinamento che a lungo termine tende a zero
curva di Kuznets in perfetto accordo con l’ipotesi per cui,
quando il PIL pro capite supera un certo valore, salvaguardia
ambientale e aumento PIL vanno nella stessa direzione.
Facendo variare l’inquinamento iniziale ed il capitale iniziale (punto
iniziale delle traiettorie), si vede che l’andamento di queste funzioni
rimane identico anche se ovviamente più è alto il valore iniziale
dell’inquinamento più tempo sarà necessario affinché si avvicini a
zero. Quindi
Investendo il 20% di capitale in Ricerca&Sviluppo si ha
crescita economica e salvaguardia dell’ambiente
Meno del 20% non è sufficiente!
A. Cutrı̀
11-01-2016, 13-01-2016
Metodi Matematici per l’ingegneria–Ing. Gestionale