Il modello HYPACT (the RAMS HYbrid Particle and Concentration Transport Model)

Il modello HYPACT (the RAMS Hybrid Particle and Concentration Transport Model)

Il modello HYPACT viene sviluppato alla Colorado State University.

 

1. Introduzione

Il modello HYPACT, disponibile presso il CRATI S.c.r.l. nella sua versione 1.0 è un codice sviluppato per simulare il moto di traccianti ed inquinanti sotto l’influenza del moto atmosferico, includendo in questo la turbolenza.
Sebbene il modello RAMS sia in grado di svolgere funzioni simili, il vantaggio fondamentale del modello HYPACT, risiede nella sua natura ibrida di trasporto.
Infatti un trattamento euleriano completo della diffusione e dispersione di traccianti, così come quello del RAMS, dovrebbe rappresentare necessariamente la sorgente con un volume non inferiore al volume di griglia e diffonderebbe rapidamente ai punti di griglia vicini la specie considerata. Uno schema lagrangiano, al contrario, consente di mantenere concentrata la piuma emessa da una data sorgente fin quando le leggi della fisica non impongono alla piuma di mescolarsi e di disperdersi.
D’altro canto uno schema di trasporto lagrangiano non è in grado di descrivere in bene il comportamento della piuma a grande distanza dalla sorgente quando quest’ultima diviene bene mescolata e dispersa. In questo caso, infatti, il numero di particelle da dover seguire per una corretta descrizione della piuma sarebbe eccessivamente dispendioso in termini di tempo computazionale e, di conseguenza, si preferisce utilizzare un approccio euleriano. La struttura ibrida del modello HYPACT consente di utilizzare un approccio lagrangiano in vicinanza della sorgente che emette la piuma ad uno schema euleriano a distanze “sufficientemente grandi” dalla sorgente.
Sebbene il modello HYPACT può ricevere informazioni da vari tipi di dati (sia da osservazioni che da campi ipotetici) il suo utilizzo naturale è in unione al modello RAMS. Nelle prossime versioni il modello HYPACT dovrebbe diventare una componente del modello RAMS e l’esecuzione del modello RAMS+HYPACT consentirà in una sola simulazione di ottenere il risultato del modello a mesoscala RAMS e del modello di trasporto ibrido HYPACT. Al momento, tuttavia, il modello HYPACT è separato dal modello RAMS e riceve da esso lo stato dinamico e termodinamico dell’atmosfera. Ciò significa che al momento eventuali processi di feedback, come ad esempio l’influenza sulla radiazione di una determinata piuma di fumo, non vengono considerati. Una volta che una certa sorgente emette in atmosfera una data specie possono intervenire su di essa eventuali termini di sorgente o pozzo dovuti a reazioni chimiche, evaporazione e condensazione, precipitazione, decadimento radiattivo etc.
I campi atmosferici vengono determinati alle varie posizioni di ogni particella rilasciata da medie spazio temporali eseguite sull’output del modello RAMS. Per ottenere un’efficienza massima nella interpolazione ed utilizzare la più alta risoluzione ottenibile dal modello RAMS, il modello HYPACT utilizza una griglia speciale chiamata “griglia – h”. La dispersione degli inquinanti nella parte lagrangiana del modello HYPACT è simulata seguendo un gran numero di particelle e riscaldando la loro posizione dalle relazioni:

(1)

(2)      

 

(3) 

dove u, v, w sono le componenti del vento direttamente risolto dal modello a mesoscala RAMS; u', v', w' sono i termini di turbolenza. wp è una eventuale velocità verticale dovuta a particolari “forcings”. Ogni particella rappresenta un elemento di massa di inquinante, wp che può variare nel tempo a causa di:

  • processi di deposizione

  • decadimento chimico o radiologico

  • evaporazione (se si considera la presenza di gocce)  

 

2. Dispersione di traccianti (con un diametro < 1 µm)

La parte non turbolenta della velocità, come detto, viene ottenuta tramite interpolazione spazio temporale delle velocità ottenute dal modello RAMS, mentre le componenti turbolente della velocità sono derivate da uno schema a catena, al primo ordine, derivato da Markov sotto le seguenti ipotesi:

  • le particelle hanno la stessa distribuzione di velocità dell’aria;

  • la turbolenza è omogenea nel piano orizzontale;

  • la turbolenza è gaussiana (ovvero le componenti turbolente del vento hanno un distribuzione normale con media nulla e una varianza nelle tre componenti delle velocità del vento);

  • le particelle si muovono indipendentemente.

La catena di Markov può essere descritta nel modo seguente:

(4)

 

(5)      

 

(6) 

  

dove ru, rv , rw sono delle componenti gaussiane non correlate, Ru , Rv , Rw sono termini di autocorrelazione lagrangiana che, nel modello HYPACT, hanno la forma semplice esponenziale.

(7)

(8)      

 

(9) 

dove su2 ,sv2 , sw2 sono le varianze delle componenti della velocità, Dt è il timestep, TLu , TLv , TLw sono i tempi di scala lagrangiani.
Le fluttuazioni ru , rv , rw sono generate ad ogni timestep dalle:

(10)

                                       

dove hu , hv, hw sono numeri casuali generati da una distribuzione gaussiana.
Le variabili meteorologiche necessarie come input per il modello HYPACT sono:

  • u, v, w  (componenti non turbolente);

  • su2 ,sv2 , sw2 (varianze del vento);

  • TLu , TLv , TLw (tempi scala lagrangiani)

Come detto il modello HYPACT ha bisogno come campi di input delle componenti non turbolente della velocità, dei tempi di scala lagrangiani e delle varianze delle velocità.
Le componenti non turbolente delle velocità vengono ottenute nel modello RAMS con semplice interpolazione spazio-temporale sulla griglia - h, mentre le varianze delle velocità ed i tempi di scala lagrangiani vengono ottenuti come descritto nel paragrafo successivo.

 

3. Parametrizzazione della turbolenza

Le variabili turbolente necessarie al modello HYPACT sono state derivate usando uno schema semplificato della tecnica elaborata da Mellor e Yamada (1974,1982).
E’ possibile applicare i seguenti schemi:

  •  Uno schema a chiusura 2.5 modificato nel caso di turbolenza crescente (Helfand e Labraga, 1988) che può essere applicato se si conoscono il campo di vento, la temperatura potenziale e l’energia cinetica turbolenta;

  •  Uno schema al secondo ordine se si hanno informazioni solo sulla temperatura e sul campo di vento. Questo schema assume un bilancio esatto tra la produzione e la dissipazione dell’energia cinetica di turbolenza.

Il vento e la temperatura entrano nei calcoli in forma di gradienti verticali adimensionali:

(11)

                                                      

La lunghezza di scala della turbolenza l è assunta da Mellor e Yamada (1982). Le varianze delle velocità ed i tempi di scala lagrangiani sono dati da:

(12)

                         

(13)

  

 

     

dove la costante cz è assunta pari a 2 (Hama, 1982).
La differenza tra lo schema a chiusura 2 e lo schema a chiusura 2.5 risiede nel calcolo dei parametri
S
m ed Shnella (12). Una volta ottenuti questi parametri si possono calcolare le varianze delle velocità dalla (12) ed i tempi scala lagrangiani dalla (13). Inoltre nel caso dello schema a chiusura 2 il calcolo dei precedenti parametri viene fatto in modo iterativo mentre per lo schema 2.5 tale calcolo è diretto.

4. Trasformazioni chimiche e decadimento radiologico

Decadimenti radiologici e chimici possono essere inclusi facilmente nel modello HYPACT deducendo la massa della specie considerata ad ogni timestep in base alla relazione:

(14)

                     

dove Dt è lo step temporale e a=(ln 2)/ t1/2  è la velocità della trasformazione .t1/2 rappresenta il tempo di dimezzamento.
Anche le trasformazioni chimiche lineari possono essere facilmente introdotte nel modello HYPACT. Ad esempio, se si considerano le specie con concentrazione
c1 e c2 , rispettivamente, che sono legate tra loro da semplici relazioni quali, ad esempio:

(15)

                    

si può facilmente risalire alla variazione della massa delle specie chimiche nel tempo dalla soluzione analitica dalle equazioni precedenti:

(16)

  

                  

Un discorso completamente diverso è invece legato alle reazioni chimiche non lineari come sottolineato in Zanetti.( Zanetti, 1992).
I coefficienti
a1, a2 possono dipendere da variabili meteorologiche quali temperatura, umidità e radiazione.
La trattazione delle reazioni chimiche non lineari, infatti, è estremamente complicata, se non impossibile, nei modelli lagrangiani.

 

5. Dispersione di particelle pesanti

La dispersione delle particelle con diametro inferiore ad un millimetro viene trattato secondo quanto illustrato nel paragrafo 2.
Nel caso di particelle con diametri superiori bisogna aggiungere gli effetti della forza gravitazionale, introducendo, nell’equazione (3) una velocità di caduta
vt . Questa velocità vale:

(17)

                    

dove Tp  è il tempo di rilassamento dato da:

(18)

                    

doverp e r sono le densità, rispettivamente, della particella e dell’aria, g è l’accelerazione di gravità, nè la viscosità molecolare cinematica. Il fattore f(Re) è introdotto per tener conto delle forze di attrito agenti sulla particella ed è funzione del numero di Reynolds:

(19)

                      

ed 

A causa delle dimensioni delle particelle in questione bisogna tener conto di altri due effetti che caratterizzano la risposta di queste particelle alla turbolenza:

  • Inerzia . Una particella con una data inerzia non risponde a tutti gli “eddies” ma solo a quelli che hanno un periodo molto maggiore del tempo di rilassamento Tp. L’effetto di inerzia può essere trascurato per particelle con diametro fino a 500mm.

  • Effetto di attraversamento degli eddies. Una particella con una data velocità di caduta non segue nel suo moto gli eddies ma continua ad attraversarli nel suo moto. Questo effetto è importante anche per particelle piccole.

L’effetto di inerzia si traduce in una diminuzione della varianza della particellaspi2 e nella crescita dell’atmosfera nel loro moto. L’effetto dell’attraversamento degli eddies si traduce in una diminuzione dell’autocorrelazione Rpi (Dt).
Di conseguenza le particelle più pesanti possono disperdersi più velocemente degli elementi di fluido se l’inerzia controlla la dispersione.
D’altro canto se l’inerzia è trascurabile , la dispersione delle particelle è inferiore a quello di un tracciante che segue il moto del fluido. Uno studio approfondito della dispersione dei traccianti pesanti si trova in Wang e Stock (1993).

 

6. Deposizione secca

L’interazione di una particella con la superficie è parametrizzata secondo Boughton et al. (1987). Sopra una certa quota h, la probabilità di deposizione è data da Monin (1959). Questa probabilità di deposizione vale:

(20)

                      

dove vd è la velocità di deposizione,f(x) è la distribuzione di Gauss standard, K è la diffusività eddy e z è l’altezza della particella al disopra della superficie.
Il calcolo della deposizione secca viene suddiviso in due casi:

  • deposizione di  particelle grandi. Le particelle che hanno z<h vengono assorbite se un numero casuale nell’intervallo chiuso [0, 1] è minore di P(z,Dt)  altrimenti la particella viene riflessa.

  •  Deposizione di particelle piccole. In questo caso si assume che la particella venga perfettamente riflessa perdendo una frazione della sua massa pari a P(z,Dt). Questa riduzione è applicata , ad esempio, nella simulazione con la chimica lineare dove la distribuzione P(z,Dt) viene calcolata ed applicata ad ogni specie chimica differente.

La velocità di deposizione vd viene calcolata con il metodo del “resistence model” (Walck et al., 1986).

 

7. Fenomeni di galleggiamento

Vari modelli sono stati suggeriti per trattare il galleggiamento delle piume emesse da ciminiere (Zanetti e Al Madami, 1984; Cofan, 1985; Gaffen et al., 1987; Van Dop, 1982) in cui una equazione stocastica descrive i fenomeni di galleggiamento della particella. Tuttavia nel modello HYPACT si preferisce adottare un approccio più semplice basato sulla derivazione di una velocità verticale addizionale che descrive, appunto, gli effetti di buoyancy. Due sono gli approcci utilizzabili nel modello HYPACT:

  • Le particelle vengono rilasciate ad una altezza effettiva zeff ;

  • Le particelle sono rilasciate al top di una ciminiera con una velocità addizionale wp che descrive il galleggiamento della particella.

La quota zeff o la velocità wp vengono ricalcolate ad ogni passo temporale quando si hanno nuovi campi meteorologici.  

 

8. Il calcolo delle concentrazioni

La massa della i-esima particella al tempo di rilascio è data da:

dove il valore di Q rappresenta l’emissione per unità di tempo nel caso di rilascio continuo ed M è la massa totale di inquinante rilasciato nell’unità di tempo. Ad esempio se rp=1/ mp per rilasciare una singola particella occorrono mp  timesteps. La sorgente da cui vengono rilasciate le particelle deve essere specificata sia nella sua posizione, all’interno del dominio considerato, sia nelle sue dimensioni. Ciò consente di poter studiare problemi relativi al rilascio da sorgente puntuale, da sorgente lineare, da sorgente areolare e da sorgente di volume.
La concentrazione
C dell’inquinante, ad un certo istante di tempo e ad una data posizione, può essere determinata contando le particelle che si trovano in un certo volume considerato
Dxs, Dys, Dzs centrato nella posizione (x,y,z):

(21)

                      

dove:

              

dove xi, yi,zi  ed mpi sono rispettivamente, la posizione e la massa della particella i-esima al tempo t.
Per poter ottenere una statistica significativa secondo uno schema appena illustrato, occorre, generalmente, rilasciare alcune migliaia di particelle. Questo rende poco efficiente il modello da un punto di vista computazionale.
Di conseguenza sono stati fatti diversi studi per migliorare questo aspetto. In particolare l’efficienza di calcolo del modello è stata migliorata con l’utilizzo di un kernel per il calcolo della concentrazione.
La concentrazione ad un dato punto viene calcolata dalla somma dei contributi di tutte le particelle considerando anche gli effetti di riflessione del suolo:

(22)

                      

dove il kernel k soddisfa la condizione:

(23)

                    

ed:

(24)

                      

I parametri hxi , hyi , hzi sono chiamati larghezze di banda. Lo schema per il calcolo della concentrazione con l’utilizzo del kernel non dipende in modo sostanziale dalla forma funzionale del kernel, tuttavia dipende in modo sensibile dalle larghezze di banda h.
Purtroppo non si è ancora arrivati alla scrittura di forme funzionali soddisfacenti per le larghezze di banda e di conseguenza, nel modello HYPACT si utilizzano le seguenti espressioni per le h:

(25)

                                                          

Come si vede dalla (25) le larghezze di banda sono costanti per tutte le particelle. La forma del kernel utilizzata dal modello HYPACT è:

(26)

                       

con:

(27)

                      

 

Bibliografia

[1]   Mellor, G.L. and T. Yamada, 1974: A hierarchy of turbulence closure models for planetary boundary layers. J. Atmos. Sci., 31, 1791-1806.

[2]   Mellor, G.L. and T. Yamada, 1982: Development of a turbulence closure model for geophysical fluid problems. Rev. Geophys. Space Phys., 20, 851-875.

[3]   Helfant, H. M., M. J. C. Labraga, 1988: Design of non singular level 2.5 second order closure model for the prediction of atmospheric turbulence. J. Atmos. Sci., 45, 113-132.

[4]   Hanna, S. R. 1982: Application in air pollution modelling. In atmospheric turbulence and air pollution modelling, Nieuwstadt, F. T. M. and H. Van Dop, editors D. Reydel publ., Dortrecth, 275-310.

[5]   Zannetti P., 1992: Particle modelling and its application for simulating air pollution phenomena. In Environmental Modelling, Melli P., and P. Zannetti editors, Computational Mechanics Publications and Elsevier Applied Science, 211-241.

[6]   Wang, L. P. and D. E. Stock, 1993: Dispersion of heavy particles by turbulent motions. J. Atmos. Sci., 50, 1897-1913.

[7]   Bougthon, B. A:, J. M. Delaurentis and W. E. Dunn, (1987): A stochastic model of particle dispersion in the atmosphere. Bound. Lay. Meteorol., 40, 147-163.

[8]   Monin A. S., 1959: On the boundary conditions on the Earth surface for diffusive pollution. Adv. Geophys., 6, 435-436.

[9]   Walchek C.J., R. A. Brost, J. S. Chang and M. L. Wesley, (1986): SO2 sulphate and HNO3 deposition velocities using regional land use and meteorological data. Atmos. Environ., 20, 949-964.

[10]  Zannetti P. and N. Al-Madani, 1984: Simulation of trasformation, buoyancy and removal processes by lagrangian particle methods. In Proc. of the 14th international technical meeting in air pollution modelling and its application, the Wispelaere, C. Editor Plenum Press, New York, 733-744.

[11]   Cogan, J. L., 1985: Monte Carlo simulation of buoyant dispersion. Atmos. Environ., 19, 867-878.

[12]  Gaffen, D. J., C. Benocci and D. Olivari, 1987: Numerical modelling of buoyancy dominated dispersion using a lagrangian approach. Atmos. Environ., 21, 1285-1293.

[13]  Van Dop, H., 1992: Buoyant plume rise in lagrangian framework. Atmos. Environ., 26A, 1345-1346.