Vanadis – model rozprzestrzeniania się zanieczyszczeń w atmosferze

W rozwiązaniu  zastosowana jest metoda Eulera, w której układ współrzędnych i siatka elementów skończonych związane są ze źródłem zanieczyszczenia. 
Zanieczyszczenie przemieszcza się względem siatki z prędkością vx, vy, vz . Przyjęty sposób dyskretyzacji obszaru za pomocą elementów ośmiowęzłowych i trójkątnych zobrazowany jest na poniższym rysunku

Dyskretyzacja obszaru trójwymiarowego (a) za pomocą elementów ośmiowęzłowych oraz obszaru dwuwymiarowego (b) za pomocą elementów trójkątnych
 
W wyniku odpowiedniego podziału obszaru na elementy uwzględniona zostaje topografia terenu. Rozwiązanie równania dyfuzji metodą elementów skończonych w przestrzeni trójwymiarowej polega na określeniu pola stężeń w węzłach elementów.
Warunki brzegowe i parametry modelu
          Równanie dyfuzji rozwiązywane jest przy przyjęciu na części powierzchni obszaru warunku brzegowego pierwszego rodzaju  S(x,y,z)=0
Warunek ten zakłada całkowity zanik zanieczyszczenia w dużej odległości od źródła emisji. Natomiast przyjęcie warunku brzegowego różnego od zera (S(x,y,z)>0) równoznaczne jest z uwzględnieniem tła zanieczyszczenia.
Na powierzchni ziemi zanieczyszczenie jest pochłaniane; gęstość pochłanianego strumienia zanieczyszczeń jest proporcjonalna do współczynnika wnikania (prędkość suchego osiadania).
W całej objętości zanieczyszczenie ulega rozpadowi z czasem połowicznego zaniku ln2/P (współczynnik zaniku P może być funkcją stężenia zanieczyszczenia, funkcją położenia i czasu oraz warunków meteorologicznych, umożliwia to uwzględnienie reakcji chemicznych, opadów i wszystkich tych czynników, które wpływają na rozpad, bądź też powstawanie zanieczyszczeń). W wybranych elementach umieszczone jest źródło zanieczyszczenia o objętościowym natężeniu emisji Qv
Numeryczne rozwiązanie równania transportu masy napotyka na trudności związane z występowaniem członu konwekcji, który może powodować do rozwiązań oscylacyjnych. Otrzymanie zadowalających wyników może być osiągnięte poprzez zmniejszenie wymiarów siatki elementów lub przez zastosowanie metody reszt ważonych z niesymetrycznymi funkcjami wagi. W tym rozwiązaniu zastosowałem metodę Galerkina.
Obliczanie wartości elementów składających się na układ równań liniowych jest realizowane metodą Gaussa. Układ równań rozwiązany jest metodą gradientów sprzężonych.
 
 
Przykładowe obliczenia (dane z pliku konfiguranyjnego vanadis.cfg):
Q[kg/s]     K[m^2/s]     P[1/s]     V1[m/s]     V2[m/s]     V3[m/s]    ALPHA[m/s]
*********************************************************************************
0         1         2         3         4         5         6         7         8            
012345678901234567890123456789012345678901234567890123456789012345678901234567890
*********************************************************************************
+001.00E+00 +040.00E+00 +07.00E-04 +00.00E+00 -00.00E+00 +01.50E+00 +01.00E-03
 
 

 
Q[kg/s]     K[m^2/s]     P[1/s]     V1[m/s]     V2[m/s]     V3[m/s]    ALPHA[m/s]
*********************************************************************************
0         1         2         3         4         5         6         7         8            
012345678901234567890123456789012345678901234567890123456789012345678901234567890
*********************************************************************************
+001.00E+00 +040.00E+00 +07.00E-04 +00.00E+00 -00.00E+00 +00.00E+00 +01.00E-03
 
 

 Porównanie modelu trójwymiarowego z modelem Pasquilla
          W celu porównania modeli przeprowadziłem obliczenia rozprzestrzeniania się zanieczyszczeń wokół źródła emisji.
Rozwiązanie Pasquilla opracowałem w oparciu o „Wytyczne obliczania stanu zanieczyszczenia powietrza atmosferycznego Ministerstwo Ochrony Środowiska, Zasobów Naturalnych i Leśnictwa, 1993” autorstwa J. Iwanka.
 
W rozwiązaniu Pasquilla przyjęte są następujące dane:
·     natężenie emisji                                          1 g/s
·     położenie źródła emisji                                 X=0, Y=0, H=60 m
·     typ równowagi atmosfery                             4
·     średnia prędkość wiatru                               1.55 m/s
·     czas połowicznego zaniku zanieczyszczenia 1000 s
·     szorstkość terenu                                        1.5 m

Rozkład zanieczyszczeń otrzymany modelem Pasquilla w płaszczyźnie XY dla Z=0 (a) oraz w płaszczyźnie YZ dla X=900m (b); vy=0
 

Rozkład zanieczyszczeń otrzymany metodą elementów skończonych w płaszczyźnie XY dla Z=0 (a) oraz w płaszczyźnie YZ dla X=900m (b); vy=0