2
∞ algoritmi numerici sofisticati, allo scopo di linearizzare e risolvere numerosi
sistemi di equazioni; possibilità di manipolazione di input e parametri non uniformi
ed eterogenei (a seconda che i dati lo permettano);
∞ regolazione di svariate condizioni al contorno del dominio studiato;
∞ la possibilità di calcolare come si ripartiscono i flussi tra scorrimento superficiale
ed infiltrazione.
In particolare questa dissertazione si soffermerà sullo studio di un modello idrologico
distribuito, efficiente su applicazioni a piccole scale, adatto a rappresentare in un
modo coerente con la realtà tutte le possibili dinamiche ed i processi di scambio idrico
tra suolo e sottosuolo.
Oltre alla descrizione teorica del modello e del suo funzionamento, si presenteranno i
principali parametri che si possono impostare al suo interno, al fine di sottolinearne la
complessità e la versatilità. Infine si cercherà di mostrare come il modello si applica,
simulando prima una situazione semplificata come caso test, ed in seguito realizzando
una calibrazione dell’idrogramma di portata simulata per un bacino canadese di 700
km
2
, di cui si è in possesso di alcuni dati reali, tra cui i valori di portata osservata alla
sezione di chiusura.
Il motivo per cui l’applicazione realizzata tratta un bacino canadese, deriva dal fatto
che questa tesi sperimentale è stata svolta all’interno di un tirocinio effettuato presso
il centro di ricerca canadese “Institut National de la Recherche Scientifique” (INRS –
ETE) della città di Quebec, nel dipartimento Acqua, Ambiente e Territorio.
3
2.
IL
MODELLO
IDROLOGICO
DISTRIBUITO
CATHY
2.1 INTRODUZIONE
CATHY (CATchment HYdrology) è un modello numerico che realizza simulazioni
accoppiate tra scorrimento superficiale e flusso idrico sotterraneo, il tutto basato
fisicamente su scale di bacino distribuite. Il tutto si basa su di un modello matematico
composto dall’accoppiamento dell’equazione di Richards per flussi tridimensionali in
mezzo poroso parzialmente saturo (modellizzato utilizzando un metodo agli elementi
finiti), unite ad una soluzione approssimata dell’equazione monodimensionale
dell’onda diffusiva per le dinamiche idrologiche superficiali.
Come tutti i modelli idrologici distribuiti, anche CATHY necessita di numerosi dati
misurati relativi al bacino analizzato, per poterli inserire nei moduli di analisi. Tra gi
input chiave troviamo dati digitalizzati del terreno, funzioni dell’idraulica del suolo,
coefficienti di scabrezza e della geometria dei ruscelli di versante e dei canali di
pianura, proprietà dell’acquifero e misure della precipitazione e dell’evaporazione
potenziale.
Il bacino idrografico viene considerato come un sistema soggetto ad ingresso di
pioggia effettiva (al netto dell’evapotraspirazione) variabile nel tempo e nello spazio,
la cui uscita (o risposta) è in parte rappresentata dall’andamento nel tempo del flusso
delle acque superficiali defluenti attraverso una o più sezioni del bacino. Oltre a ciò, il
modello calcola diverse variabili di stato distribuite nello spazio e nel tempo (grado di
saturazione, flussi sotterranei, flussi di interscambio tra superficie e sottosuolo,
variazione del livello di falda).
CATHY in particolare si compone di un modulo di analisi sotterranea (descritto
dall’equazione di Richards), accoppiato ad un modulo di analisi superficiale (descritto
dall’equazione dell’onda diffusiva), integrato da un processo di analisi dei dati
territoriali digitalizzati al fine di realizzare una precisa delineazione del reticolo di
drenaggio; la determinazione della rete idrografica risulta di fondamentale
importanza, in quanto definisce dove è presente il deflusso canalizzato e dove si
verifica quello di versante, due tipologie di movimento idrico completamente diverse
e regolate da processi fisici differenti. Nel modello inoltre vengono presi in
4
considerazione gli effetti di accumulo e ritardo dettati dalla presenza di laghi o altre
depressioni della topografia, che influenzano le dinamiche dei flussi idrici.
L’accoppiamento tra i moduli di analisi superficiale e del sottosuolo è basato su una
procedura di variazione dinamica (“switching”) delle condizioni al contorno, in cui
viene considerata l’altezza dell’acqua sulla superficie per suddividere i flussi
potenziali atmosferici in ingresso tra flussi effettivi di scorrimento superficiale ed
immagazzinamento idrico. Le condizioni iniziali sono generate a partire da
informazioni di conduttività superficiale e di altezza della lama d’acqua sul suolo,
mentre gli intervalli di tempo per il calcolo della simulazione sono variabili.
2.2 DESCRIZIONE DEL MODELLO
Il modello CATHY descrive i processi idrologici tramite un modello matematico
rappresentato da un accoppiamento tra una soluzione agli elementi finiti della
equazione di Richards per flussi in mezzo parzialmente saturo, ed una soluzione alle
differenze finite dell’equazione dell’onda diffusiva, che rappresenta la propagazione
dello scorrimento superficiale attraverso una connessione tra corsi d’acqua di versante
e di canale (classificati come tali utilizzando un concetto di topografia del terreno e di
geometria idraulica).
Di seguito si riporta il modello matematico, descritto da un sistema di equazioni
differenziali alle derivate parziali (entrambe basate sul principio di conservazione
della massa e del momento):
€
σ S
w
( )
∂ψ
∂t
= ∇ ⋅ K
s
K
rw
S
w
( )
∇ψ + η
z
( )
[ ]
+ q
ss
h
( )
∂Q
∂t
+ c
k
∂Q
∂s
= D
h
∂
2
Q
∂s
2
+ c
k
q
s
h,ψ
( )
5
Dove per la prima equazione:
∞
€
σ S
w
( )
= S ψ
( )
=
∂θ
∂ψ
è il termine per indicare l’immagazzinamento idrico
(“storativity”)
€
L
−1
[ ]
;
∞
€
S
w
=
θ
θ
s
è il grado di saturazione;
∞
€
θ è il contenuto volumetrico di acqua nel suolo (umidità del suolo);
∞
€
θ
s
è l’umidità del suolo a saturazione, che coincide con la sua porosità
€
φ ;
∞
€
ψ = p /γ = γ h /γ = h è un valore di pressione (potenziale capillare)
€
L
[ ]
;
∞
€
t e il tempo
€
T
[ ]
;
∞
€
∇ è l’operatore di gradiente;
∞
€
K
s
è il tensore della conduttività idraulica in condizioni di saturazione
€
L / T
[ ]
;
∞
€
K
rw
S
w
( )
= K
r
ψ
( )
è la funzione della conduttività idraulica relativa (adimensionale);
∞
€
η
z
= 0,0,1
( )
T
;
∞
€
z è la coordinata verticale (positiva verso l’alto)
€
L
[ ]
;
∞
€
q
ss
h
( )
rappresenta i termini distribuiti di sorgente (positivi) o di assorbimento
(negativi) dalla superficie al sottosuolo
€
L
3
/ L
3
T
[ ]
.
Mentre per la seconda equazione:
∞ per rappresentare i singoli collegamenti tra i versanti ed i canali viene utilizzato un
sistema di coordinate curvilinee,
€
s
€
L
[ ]
;
∞
€
Q è il valore di portata lungo il canale/versante
€
L
3
/ T
[ ]
;
∞
€
c
k
è la celerità dell’onda cinematica
€
L / T
[ ]
;
∞
€
D
h
è la diffusività idraulica
€
L
2
/ T
[ ]
;
∞
€
q
s
h ,ψ
( )
= i − f
i
+ f
e
− e
p
( )
A /L è il flusso di sorgente (positivo) o di assorbimento
(negativo) tra il sottosuolo e il terreno superficiale
€
L
3
/ LT
[ ]
, in cui:
⇒
€
i è l’intensità della pioggia
€
L /T
[ ]
;
⇒
€
f
i
è il flusso effettivo di infiltrazione
€
L /T
[ ]
;
⇒
€
f
e
è il flusso effettivo di exfiltrazione
€
L /T
[ ]
;
⇒
€
e
p
la domanda atmosferica di evaporazione
€
L /T
[ ]
;
6
⇒
€
A è l’area di terreno superficiale su cui sussiste l’elemento idrico analizzato
€
L
2
[ ]
;
⇒
€
L è la lunghezza del singolo elemento analizzato della rete idrica
€
L
[ ]
.
Figura 1 - Schema rappresentativo della dinamica di scambio dei flussi idrici in CATHY
[Fonte: Camporese et al. (2008)]
Questo sistema di equazioni deve essere risolto simultaneamente per il vettore di
incognite
€
Q,ψ
( )
o
€
h ,ψ
( )
. Le non linearità crescono nelle curve caratteristiche
dell’equazione di Richards
€
S
w
ψ
( )
e
€
K
rw
S
w
( )
, a causa della dipendenza dalle variabili
non lineari di
€
q
ss
per
€
h e
€
q
s
per
€
ψ . Questo necessita una risoluzione con un metodo di
linearizzazione utilizzando l’iterazione di Picard; a questa si aggiunge una divisione
automatica degli intervalli di tempo di analisi delle equazioni, con un controllo
dell’errore integrato utile ad aumentare la durata media degli intervalli temporali,
incrementando di conseguenza l’efficienza computazionale totale del modello.
Con questo modello è così possibile gestire e simulare automaticamente afflussi e
deflussi idrici, che possono essere alternativamente in funzione sia delle condizioni
atmosferiche che delle caratteristiche del suolo; inoltre si possono rappresentare le
due tipologie di formazione di scorrimento superficiale, cioè per eccesso di
saturazione (secondo Dunne) o per eccesso di infiltrazione (secondo Horton).
7
2.2.1 Modulo Superficiale
La modellazione dei fenomeni idrologici superficiali in CATHY avviene attraverso la
definizione di due processi legati a flussi di canale o di versante, a cui si applica una
differente parametrizzazione. Una volta distinti i flussi precedenti, si crea un bacino di
drenaggio che rappresenti la rete idrica suddivisa tra i corsi d’acqua di versante
(rivoletti) ed i canali, in cui siano considerati anche gli effetti di immagazzinamento e
ritardo nei fiumi e nei laghi (a cui si aggiungono gli effetti di
infiltrazione/exfiltrazione dal sottosuolo).
I rivoletti ed i canali sono quindi rappresentati entrambi attraverso una stessa
equazione monodimensionale di diffusione-advezione, in cui si utilizzano differenti
valori di parametri, come descritto da Orlandini e Rosso (1996 e 1998), finalizzati a
distinguere le geometrie idrauliche dei due differenti regimi di corrente.
Per descrivere la composizione e la topografia del sistema drenante, precedentemente
al lancio del modello di analisi, viene estratta automaticamente una rete di drenaggio a
partire dai modelli digitalizzati delle quote del terreno DEM (Digital Elevation
Model). In questa procedura ad ogni elemento si associano una lunghezza ed una
pendenza locale del terreno (che dipenderanno dalla cella corrispondente estratta
automaticamente dalla rete di drenaggio), oltre che una larghezza di sezione
trasversale rettangolare (la cui ampiezza varia dinamicamente con la portata in uscita,
secondo le proprietà della geometria idraulica alle diverse scale). La geometria
idraulica dei canali nella modellazione della diffusione delle onde per le dinamiche di
bacino, è sintetizzata combinando le relazioni fluviali “at-a-station” e “downstream”
di Leopold e Maddock (1953), relative alla larghezza del pelo libero ed al perimetro
bagnato.
Il modello superficiale ha il compito di calcolare l’altezza della lama d’acqua
€
h
che, a
seconda che sia maggiore o uguale a zero, determina la condizione atmosferica al
contorno (rispettivamente di Dirichlet o Neumann), per impostare successivamente
l’analisi del modulo sotterraneo.
2.2.1.1 Estrazione della rete di drenaggio
Nell’algoritmo utilizzato per delineare il reticolo di drenaggio è possibile scegliere
uno dei tre metodi proposti in letteratura da O’Callaghan e Mark (1984), da Tarboton
(1997) e Orlandini et al. (2003), conosciuti col nome di
€
D8,
€
D8 − LAD e
€
D8 − LTD.
8
Nello specifico tutti e tre gli algoritmi si caratterizzano nell’assegnare una direzione di
scorrimento del flusso che passa per il centro di ogni singola cella quadrata, verso il
centro delle otto celle adiacenti che la contornano seguendo le quattro direzioni
cardinali e le quattro diagonali (con conseguenti approssimazioni). In particolare le
più recenti metodologie
€
D8 − LAD e
€
D8 − LTD risultano più precise, poiché
introducono una suddivisione ulteriore in triangoli dello spazio circostante la cella
analizzata, i quali aiutano nel ponderare la teorica direzione del flusso, da scegliere tra
due celle confinanti.
Infatti per quanto riguarda il metodo
€
D8 − LAD (Least Angular Deviation), cioè della
minor deviazione angolare, la scelta della cella di drenaggio successiva a quella
analizzata la si effettua in base ad una ponderazione dell’ampiezza dell’angolo (
€
α ) tra
la direzione teorica del flusso e quella diagonale o cardinale, che deve essere minore.
Mentre invece per il metodo
€
D8 − LTD (Least Transversal Deviation), cioè della
minor deviazione trasversale, la scelta della cella di drenaggio successiva a quella
analizzata avviene in base ad una minimizzazione della distanza lineare cumulativa
(
€
δ) tra il centro della ipotetica cella drenante e la direzione teorica del flusso.
Nella figura sottostante è possibile comprendere meglio quanto ora enunciato.
Figura 2 – Schema rappresentativo del funzionamento dei metodi
€
D8 − LTD e
€
D8 − LAD
[Fonte: Orlandini et al. (2003)]