46
Komputerowa symulacja
rozprzestrzeniania
zanieczyszczeń w atmosferze
Jacek Piekarski, Marcin Lubierski
Katedra Techniki Hydro-Szlamowej
i Utylizacji Odpadów
Politechnika Koszali
ń
ska
1. Wstęp
Atmosfera otaczająca powierzchnię Ziemi dzięki swojej naturze posia-
da zdolność zmniejszania i kumulowania stężenia oraz przemieszczania emito-
wanych zanieczyszczeń, które następnie ulegają licznym przemianom fizyko-
chemicznym. Bezpośredni wpływ na te przemiany ma zarówno cyrkulacja
obejmująca kontynent, a nawet całą planetę, jak i nakładające się nań różno-
rodne systemy prądów powietrznych mniejszej skali. Znajomość procesów
zachodzących w atmosferze, warunków meteorologicznych panujących na da-
nym obszarze jest niezbędna dla prawidłowego przewidywania skutków emisji
zarówno ze źródeł istniejących jak i projektowanych [6,10].
Na potrzeby prognozowania rozprzestrzeniania stężeń zanieczyszczeń
opracowano szereg modeli matematycznych analityczno-empirycznych. Modele
wykorzystywane są przy tworzeniu systemów alarmowych oraz do specyfikacji
obszarów zagrożeń na wypadek wystąpienia dużych emisji losowych. Można je
stosować do symulacji występujących stężeń w atmosferze, co jest wskazane ze
względu na wysokie koszty odpowiednich pomiarów. Dlatego też znając para-
Jacek Piekarski, Marcin Lubierski
VI Ogólnopolska Konferencja Naukowa
666
metry emitora, warunki meteorologiczne i terenowe, z wystarczającą dokładno-
ś
cią można określić stężenia zanieczyszczeń w emisji.
Można wyróżnić trzy podstawowe kategorie modeli: Eulera, Lagrange'a
oraz Gaussa. W modelu Eulera rozpraszanie się zanieczyszczeń jest opisywane
względem nieruchomego układu związanego z ziemią. Zmienne w przestrzeni
i czasie, pole prędkości przedstawione jest w każdym punkcie badanego obsza-
ru, w ustalonym na powierzchni ziemi układzie współrzędnych kartezjańskich
lub sferycznych. Eulerowskie pole prędkości obrazują izogony, izotachy albo
linie prądu. Gdy jest ono stacjonarne, to linie prądu pokrywają się z trajekto-
riami cząstek. W takim polu następstwo poruszających się cząsteczek powietrza
jest widziane przez nieruchomego obserwatora jako przenoszone przez wiatr.
W modelu Lagrange'a zmiany stężenia zanieczyszczenia są opisywane
w układzie związanym z przemieszczającymi się masami powietrza. Obserwa-
tor poruszający się wraz z powietrzem utrzymuje kontakt z tymi samymi cząst-
kami powietrza w zadanym okresie czasu. W modelu gaussowskim stężenie
zanieczyszczenia wzdłuż osi smugi ma statystyczny rozkład Gaussa. Formuła ta
może być wyprowadzona przy zastosowaniu metody Eulera lub Langrage`a.
Jednak modele gaussowskie w swej najbardziej uproszczonej postaci winne być
postrzegane jako specjalne uproszczenie metody Langrange`a.
Obok trzech podstawowych modeli (Eulera, Lagrange`a i Gaussa), ist-
nieją również inne określające rozprzestrzenianie się zanieczyszczeń w atmos-
ferze [2,3].
Proste modele statystyczne stosowane głównie w praktyce inżynierskiej
wykorzystują modele regresyjne, począwszy od prostych tablic rozdzielczych
i regresji wielokrotnej do analizy ciągów chronologicznych. Po zainstalowaniu
sieci pomiarowej monitoringu, pojawiły się w Polsce proste modele progno-
styczne dyspersji zanieczyszczeń, stosowane w skali lokalnej.
Opracowany przez polskich matematyków model stochastyczny progno-
zuje pole stężeń tlenków azotu, ulegających w atmosferze złożonym przemianom.
Ponieważ charakterystyczną cechą tlenków azotu jest losowe powstawanie oraz
zanikanie, dlatego trajektorie cząstek NO
x
są nieciągłe. Niniejsze zjawiska opisy-
wane są w dwóch przestrzeniach probabilistycznych. Jedna dotyczy procesów
mikroskopowych, natomiast druga procesów makroskopowych.
Lagrangeowskie modele rozprzestrzeniania cząstek są obliczeniowo
kosztowne, dlatego w przypadku większych obszarów, wykorzystuje się modele
hybrydowe lagrangeowsko–eulerowskie (modele dyfuzyjno-stochasty-czne).
W otoczeniu punktowego źródła emisji obliczenia przeprowadzane są modelem
LPD (lagrangeowski model cząstek). Po dostatecznie długim czasie przebiegu,
gdy smuga ma duże rozmiary w porównaniu z przestrzennym krokiem całko-
Komputerowa symulacja rozprzestrzeniania zanieczyszcze
ń
w atmosferze
VI Ogólnopolska Konferencja Naukowa
667
wania, cząstki włączają się do objętościowego pola emisji modelu eulerowskie-
go i zanikają.
Na podstawie Rozporządzenia Ministra Środowiska z dnia 5 grudnia
2002 r., do określania stanu zanieczyszczeń powietrza dla źródeł istniejących
stosuje się model gaussowski (formułę Pasquilla)[8].
Pasquill przystosował uproszczone równanie dyfuzji zanieczyszczenia
gazowego w poruszającym się ośrodku gazowym, do zagadnień rozprzestrze-
niania się zanieczyszczeń wyrzucanych do atmosfery z kominów. Równanie
dyfuzji w postaci Pasquilla ma postać [4,6,7,9]:
2
2
2
2
2
2
z
S
2
Sz
dt
d
y
S
2
Sy
dt
d
t
S
∂
∂
⋅
⋅
+
∂
∂
⋅
⋅
=
∂
∂
(1)
gdzie:
S – stężenie zanieczyszczenia w punkcie recepcyjnym,
t – czas wędrówki zanieczyszczenia od chwili wyrzutu do osiągnięcia
punktu recepcyjnego,
Sy – współczynnik dyfuzji Pasquilla wzdłuż osi 0Y w układzie współ-
rzędnych t, y, z,
Sz – współczynnik dyfuzji Pasquilla wzdłuż osi 0Z w układzie współ-
rzędnych t, y, z.
Równanie opisuje rozprzestrzenianie się zanieczyszczeń dla dowolnego
punktu emisji, ustalonego w czasie. Może to być źródło punktowe, liniowe lub
powierzchniowe – zależy to od przyjętych warunków brzegowych:
ź
ródło emisji o masowym natężeniu E jest umieszczone w punkcie (źródło
punktowe) o współrzędnych t = 0, y = 0 i z = H,
na powierzchni podłoża nie zachodzi pochłanianie zanieczyszczeń, tj. dla
z = 0 strumień zanieczyszczeń jest zerowy (ilość zanieczyszczeń docierają-
cych do podłoża jest równa ilości zanieczyszczenia od podłoża),
stężenie zanieczyszczeń w odległości nieskończenie dużej od źródła (tzn. t,
y oraz z dążą do nieskończoności) ma wartość zerową.
Jacek Piekarski, Marcin Lubierski
VI Ogólnopolska Konferencja Naukowa
668
Równanie spełniające warunki brzegowe ma postać:
(
)
(
)
Sz
Sy
Usr
0063
,
0
)
e
e
(
e
Eg
S
2
2
2
2
2
2
Sz
2
Hp
z
Sz
2
Hp
z
Sy
2
y
⋅
⋅
⋅
+
⋅
⋅
=
+
−
−
−
−
(2)
gdzie:
S – stężenie substancji gazowej w punkcie (x, y, z), [µg/m
3
],
Eg – maksymalna emisja substancji gazowej z źródła, [mg/s],
z – wysokość, dla której oblicza się stężenie substancji zanieczyszczają-
cej w powietrzu, [m],
Hp – efektywna wysokość emitora, [m],
Usr – średnia prędkość wiatru w warstwie powietrza od z = h do z = Hp,
[m/s],
h – geometryczna wysokość emitora, [m],
Sy – współczynnik poziomej dyfuzji atmosferycznej, [m],
Sz – współczynnik pionowej dyfuzji atmosferycznej, [m].
Wykorzystując zjawisko, że główny ruch zanieczyszczenia z prędkością
wiatru Usr odbywa się wzdłuż osi x i zawsze występuje określona relacja po-
między tymi wielkościami, w wzorze (2) w miejsce zmiennej t podstawia się
współrzędną x.
Smuga gazów wylotowych wydobywa się z punktu pozornej emisji A
znajdującego się na wysokości Hp [m] – rysunek 1. Efektywna wysokość emi-
tora Hp stanowi sumę geometrycznej wysokości komina h [m] i tzw. wysokości
wyniesienia Dh [m]. Wartość wyniesienia Dh zależna jest od prędkości wylo-
towej gazów v [m/s], emisji ciepła Q [kJ/s] i prędkości wiatru Uh [m/s] na wy-
sokości wylotu z emitera. W zależności od ilości emitowanego ciepła Q wyróż-
niamy dwa rodzaje wyniesienia: termiczne oraz dynamiczne. Wartość wynie-
sienia termicznego obliczana jest na podstawie wzoru CONCAWE, który został
opracowany dla kominów energetycznych. Zgodnie z Rozporządzeniem Mini-
stra Środowiska z dnia 5 grudnia 2002 r. [8], granicę jego zastosowania okre-
ś
lono dla wartości emisji ciepła powyżej 24000 [kJ/s]. Formuła Hollanda (obli-
cza wysokość wyniesienia dynamicznego gazów) ma zastosowanie, gdy ilość
emisji ciepła jest mniejsza od 16000 [kJ/s]. Uwzględnia ona różnicę pomiędzy
prędkością wiatru na wysokości emitora Uh a prędkością gazów odlotowych v.
W przypadku, gdy wielkość emisji ciepła mieści się w przedziale od 16000 [kJ/s]
do 24000 [kJ/s] stosowane są obie formuły, pomniejszone o różnicę pomiędzy
Komputerowa symulacja rozprzestrzeniania zanieczyszcze
ń
w atmosferze
VI Ogólnopolska Konferencja Naukowa
669
wartością emisji ciepła Q a wartościami granicznymi. Powyższe wzory na wynie-
sienie gazów wylotowych zakładają pionowy wyrzut gazów do atmosfery. Stąd
przy wyrzucie poziomym, w wyniku zadaszenia czy wyrzutni poziomej należy
przyjąć zerową wartość wyniesienia.
Rys. 1.
Graficzna interpretacja formuły Pasquilla (opis w tekście) [6]
Rys. 1. The graphic interpretation of Pasquill formula (description in text) [6]
Wymiar poprzeczny smugi zanieczyszczeń stanowią półosie elipsy
o powierzchni, przez którą przepływa 90% masy wydalanej z punktowego źró-
dła emisji, w tym także 90% określonego zanieczyszczenia – rysunek 1. Przy
takim założeniu półosie elipsy wynoszą 2,15 Sy i 2,15 Sz. Współczynniki dyfu-
zji Sy i Sz mają wymiar długości. Są one proporcjonalne do wymiarów po-
przecznych rozszerzającej się i niezakłóconej powierzchnią ziemi smugi. Wraz
ze wzrostem odległości od źródła emisji wartości współczynników dyfuzyjnych
rosną ze względu na rozpraszanie się smugi. Istnieje kilka metod obliczania ich
wartości, najczęściej stosowaną w Polsce jest metoda Nowickiego. Charakter
zmian współczynników opisany jest przez wyrażenia [1,4,8]:
a
x
A
Sy
⋅
=
(3)
Jacek Piekarski, Marcin Lubierski
VI Ogólnopolska Konferencja Naukowa
670
b
x
B
Sz
⋅
=
(4)
gdzie:
A, B, – współczynniki empiryczne, [m],
a, b – współczynniki empiryczne, [–].
Wykładniki potęgowe a i b charakteryzują dynamikę rozpraszania smu-
gi zanieczyszczeń w atmosferze. Maleją wraz ze zmniejszaniem się intensyw-
ności turbulencji – dotyczy to szczególnie wykładnika b decydującego o kącie
pionowego rozwarcia smugi. Zależne są od stanu równowagi atmosfery, przy
czym nie są one opisane w postaci reguły matematycznej, lecz ich wartość do-
biera się na podstawie tablic [4,6].
Parametry A i B zależne są od warunków terenowych, w których pro-
wadzony był pomiar współczynników dyfuzji atmosfery. Stosunek pomiędzy
efektywną wysokością emitera a aerodynamiczną szorstkością terenu jest kryte-
rium opisującym wpływ czynników terenowych na wielkość dyfuzji. Uwzględ-
niony jest również wpływ stanu równowagi atmosfery poprzez wykładnik me-
teorologiczny m dobierany podobnie jak wykładniki a i b [4]:
)
0
Z
Hp
ln
1
m
6
(
088
,
0
A
3
,
0
−
+
⋅
=
−
(5)
)
0
Z
Hp
ln
7
,
8
(
m
38
,
0
B
3
,
1
−
⋅
=
(6)
gdzie:
m – współczynnik empiryczny, zależny od stanu atmosfery, [–],
Hp – efektywna wysokość emitora, [m],
Z0 – aerodynamiczny współczynnik szorstkości terenu, [m].
W modelu uwzględnione są (rysunek 2): prędkość wiatru na wysokości
emitora Uh [m/s], średnia prędkość wiatru w warstwie powietrza od z = h do
z = Hp Usr [m/s] oraz prędkość wiatru na wysokości anemometru Ua [m/s].
Prędkość wiatru Uh, obliczana jest w zależności od wysokości emitora h, pręd-
kości wiatru na wysokości anemometru Ua oraz stanu atmosfery (parametru
meteorologicznego m). Natomiast wartość średniej prędkości Usr zależy od
efektywnej wysokości emitora Hp, wysokości emitora h, prędkości wiatru na
wysokości anemometru Ua oraz stanu atmosfery (parametr m).
Komputerowa symulacja rozprzestrzeniania zanieczyszcze
ń
w atmosferze
VI Ogólnopolska Konferencja Naukowa
671
U
Usr
Uh
Ua
A
X
Z
1
4
H
p
h
D
h
anemometr
Rys. 2.
Rozkład prędkości wiatru wraz z wysokością w modelu Pasquilla [6]
Rys. 2. The distribution of wind velocity of along with the height in Pasquill model [6]
Aerodynamiczny współczynnik szorstkości terenu Z0 ma wymiar dłu-
gości i uwzględnia wpływ pokrycia powierzchni na intensywność rozpraszania
zanieczyszczeń w atmosferze. Elementy o zróżnicowanej wysokości mają bez-
pośredni wpływ na poziomy ruch mas powietrza, rozdzielając je na poszcze-
gólne strumienie i powodując zawirowania turbulencyjne.
Algorytm umożliwiający obliczenie stężenia zanieczyszczeń S, współczyn-
ników dyfuzji Sy i Sz, wysokości Hp pozornego punktu emisji A oraz pozostałych
niezbędnych parametrów emitora i meteorologicznych przedstawia rys. 3.
2. Opis programu
Na podstawie algorytmu przedstawionego na rysunku 3 wykonano pro-
gram komputerowy „Atmo” służący do symulacji rozprzestrzeniania się zanie-
czyszczeń w atmosferze. Program pracuje w środowisku Windows
95/98/XP/NT. Składa się z dwóch podstawowych modułów projektowych
umożliwiających wprowadzenie danych wstępnych i wyprowadzenie wyników
obliczeń oraz dwóch modułów pomocniczych służących wizualizacji przepro-
wadzonej symulacji poprzez przedstawienie zmian wartości współczynników
dyspersji oraz stężenia w zależności od zmiany wartości odległości receptora
od emitera. Aplikacja posiada opcję zapisu do pliku oraz odczytu z pliku
wprowadzonych danych początkowych. Pliki skojarzone z aplikacją mają roz-
szerzenie „*.zan”. Aplikacja generuje raport końcowy oraz daje możliwość
wydruku wprowadzonych danych początkowych i wyników obliczeń.
Jacek Piekarski, Marcin Lubierski
VI Ogólnopolska Konferencja Naukowa
672
h, d, v, T, Z0, Ua, T0, ATM, x, y, z, Eg
h > 300
Uh = 21,429·Ua
m
Uh = Ua·(0,071·h)
m
Uh < 0,5
Uh = 0,5
T
T
N
T
N
v > 0,5Uh
DhH = 0
v < Uh
Q > 16000
Dh=DhH
Q > 24000
Hollanda
CONCAWE
Hollanda + CONCAWE
Hp=h+Dh
N
N
T
T
T
N
N
Hp = h
Hp > 300
h < 300
Usr = Ua·(0,071·h)
m
Usr = 21,429·Ua
m
Usr < 0,5
Usr = 0,5
Sy = A·x
a
Sz = B·x
b
Uh
Q
0,00974
v
1,5
DhH
⋅
+
⋅
=
0,5Uh
0,5Uh
v
Uh
Q
0,00974
v
1,5
DhH
−
−
⋅
+
⋅
=
0,7
0,58
Uh
Q
1,126
DhC
⋅
=
8000
16000
-
Q
Uh
Q
1,126
8000
Q
-
24000
DhH
Dh
0,7
0,58
⋅
⋅
+
⋅
=
(
)
(
)(
)
m
m
1
m
1
14
1
m
h
-
Hp
h
-
Hp
Ua
Usr
⋅
+
⋅
=
+
+
(
)
(
)
⋅
−
+
+
−
⋅
⋅
=
+
+
m
m
1
m
1
m
300
300
Hp
m
1
h
300
14
h
-
Hp
Ua
Usr
10
0
Z
Hp
<
10
0
Z
Hp
=
1500
0
Z
Hp
>
1500
0
Z
Hp
=
−
+
⋅
=
−
0
Z
Hp
ln
1
m
6
088
,
0
A
3
,
0
−
⋅
⋅
=
0
Z
Hp
ln
7
,
8
m
38
,
0
B
3
,
1
(
)
(
)
Sz
Sy
Usr
0063
,
0
Sz
2
Hp
z
exp
Sz
2
Hp
z
exp
)
Sy
2
y
exp(
Eg
S
2
2
2
2
2
2
⋅
⋅
⋅
+
−
+
−
−
⋅
−
⋅
=
S
K
T
T
T
N
N
N
N
N
N
T
T
T
T
T0)
-
·v·(T
278,76·d
Q
2
=
m, a, b, Uh, Q, Dh, Hp, Usr, A, B, Sy, Sz, S
Rys. 3.
Algorytm umożliwiający obliczenie stężenia zanieczyszczeń S (opis w tekście)
Rys. 3. The algorithm enabling calculation of pollutants concentration S (description in
text)
Komputerowa symulacja rozprzestrzeniania zanieczyszcze
ń
w atmosferze
VI Ogólnopolska Konferencja Naukowa
673
Rys. 4.
Wprowadzanie danych
Rys. 4. Data input
Po wybraniu z menu „Zanieczyszczenia” opcji „Dane wstępne” (F5) na
ekranie monitora wyświetlone zostaje okno przedstawione na rysunku 4.
Moduł umożliwiający wprowadzenie danych projektowych podzielony
jest na cztery części charakteryzujące: parametry techniczne emitera, parametry
meteorologiczne, lokalizację receptora oraz tzw. wartości inne. Do parametrów
technicznych zalicza się: wysokość emitera h [m], średnicę wewnętrzną d [m]
wylotu przewodu emitującego substancje, prędkość v [m/s] i temperaturę T [K]
gazów na wylocie z emitera oraz maksymalną emisję Eg [mg/s] substancji ga-
zowej. Emisję substancji zanieczyszczających ze źródeł projektowanych ustala
się szacunkowo na podstawie wskaźników emisji tychże substancji, charaktery-
stycznych dla procesu technologicznego lub przez analogię do emisji ze źródła
istniejącego. W przypadku źródeł istniejących i źródeł projektowanych o przy-
jętej technologii i rozwiązaniach technicznych, emisję substancji zanieczysz-
Jacek Piekarski, Marcin Lubierski
VI Ogólnopolska Konferencja Naukowa
674
czających ustala się na podstawie wyników pomiarów, założeń technologicz-
nych oraz wielkości produkcji.
Parametry meteorologiczne to: prędkość wiatru Ua [m/s] mierzona na
wysokości anemometru, średnia temperatura powietrza T0 [K] odnosząca się
do okresu obliczeniowego oraz stan równowagi atmosfery ATM [–]. Lokaliza-
cja receptora podawana jest w formie odległości od emitera: x [m] (równoległa
z kierunkiem wiatru) i y [m] (prostopadła do kierunku wiatru) oraz wysokości,
dla której określa się stężenie substancji w terenie z [m]. Parametry inne to
współczynnik aerodynamicznej szorstkości terenu Z0 [m], przy pomocy którego
uwzględnia się wpływ pokrycia powierzchni na intensywność rozpraszania
zanieczyszczeń w atmosferze. Jego wartość wyznacza się na podstawie map
topograficznych.
Podczas edycji poszczególnych wartości parametrów obliczeniowych
program automatycznie dokonuje sprawdzenia ich poprawności, pod względem
przynależności do określonych zakresów, w których równania opisujące zmiany
stężenia zanieczyszczeń w atmosferze są słuszne. Błędy sygnalizowane są od-
powiednim komunikatem ukazującym się w oknie modalnym.
Po prawidłowym wprowadzeniu danych wstępnych należy wybrać
z menu „Zanieczyszczenia” opcję „Wyniki obliczeń” (F6) – rysunek 5, program
zacznie realizować algorytm przedstawiony na rysunku 3.
Na podstawie wartości stanu równowagi atmosfery ATM [–], która
zmienia się w granicach od 1 do 6, program dobiera parametry meteorologiczne
m [–], a [–] oraz b [–], które są używane w dalszych obliczeniach. W zależności
od wysokości h [m] emitora oblicza się prędkość wiatru Uh [m/s] na wysokości
wylotu z emitora. W przypadku, gdy obliczona wartość będzie mniejsza zosta-
nie przyjęta Uh = 0,5 [m/s]. Emisja ciepła Q [kJ/s] z emitera obliczana jest
z zależności od średnicy d [m] emitora, prędkości v [m/s] gazów na wylocie
z emitera oraz temperatur T0 [K] otoczenia i T [K] gazów na wylocie z emitera.
Następnie na podstawie wartości emisji ciepła Q [kJ/s] z emitera program do-
biera i wyświetla nazwę formuły, według której następnie oblicza wartość pa-
rametru wyniesienia Dh [m] gazów odlotowych. Dla Q < 16000 [kJ/s] zastoso-
wana zostaje formuła Hollanda, w której w zależności od stosunku prędkości
wiatru Uh [m/s] na wysokości emitera do prędkości v [m/s] gazów na wylocie z
emitera jest obliczane Dh, jako funkcja emisji ciepła Q, średnicy d emitora oraz
prędkości v gazów. Dla Q > 24000 [kJ/s] stosowana jest formuła CONCAWE,
gdzie wartość wyniesienia Dh jest obliczana na podstawie wartości emisji
Q i prędkości wiatru Uh. W przypadku, gdy Q będzie zawierać się w przedziale
wartości od 16000 do 24000 [kJ/s] zastosowane zostaną obie formuły, odpo-
wiednio pomniejszone o różnicę pomiędzy wartością emisji Q a wartościami
granicznymi. Suma wartości parametrów wysokości emitera h [m] oraz wynie-
Komputerowa symulacja rozprzestrzeniania zanieczyszcze
ń
w atmosferze
VI Ogólnopolska Konferencja Naukowa
675
sienia Dh [m] gazów odlotowych stanowi efektywną wysokość emitera Hp [m].
W zależności od efektywnej wysokości Hp emitora oraz wysokości h emitora
program oblicza średnią prędkość wiatru Usr [m/s] w warstwie, której wyso-
kość waha się od z = h do z = Hp. Jej wartość jest obliczana na podstawie efek-
tywnej wysokości emitera Hp, wysokości emitera h oraz prędkości wiatru Ua.
Program, podobnie jak dla Uh, przyjmie Usr = 0,5 [m/s], gdy średnia prędkość
wiatru będzie mniejsza. Następnie dokonywana zostaje weryfikacja ilorazu
efektywnej wysokości emitera Hp [m] oraz współczynnika aerodynamicznej
szorstkości terenu Z0 [m], który powinien zawierać się w przedziale wartości
od 10 do 1500.
Rys. 5.
Wyprowadzenie wyników symulacji komputerowej
Rys. 5. Output results of computer simulation
Jacek Piekarski, Marcin Lubierski
VI Ogólnopolska Konferencja Naukowa
676
W kolejnym etapie realizacji algorytmu (rysunek 3), program
oblicza współczynniki A [–] oraz B [–], które są zależne od współczynnika
meteorologicznego m oraz stosunku efektywnej wysokości Hp emitora
i współczynnika aerodynamicznej szorstkości terenu Z0. Na podstawie obliczo-
nych współczynników A i B, dobranych wcześniej współczynników empirycz-
nych a i b oraz odległości x [m] receptora od emitora określane zostają współ-
czynniki dyspersji poziomej Sy [m] i pionowej Sz [m]. Ostatnim parametrem,
wyznaczanym przez program jest stężenie jedno-godzinne substancji gazowej
S [µg/m
3
], zależne od wartości emisji Eg [mg/s] gazów, średniej prędkości Usr
wiatru, odległości y [m] emitora od receptora, wysokości z [m] receptora, efek-
tywnej wysokości emitora Hp oraz współczynników dyspersji Sy i Sz.
Rys. 6.
Graficzne przedstawienie zmiany współczynnika dyspersji w zależności od
odległości
Rys. 6. Graphic represaentation of dispersion coefficient change in dependence on
distance
Komputerowa symulacja rozprzestrzeniania zanieczyszcze
ń
w atmosferze
VI Ogólnopolska Konferencja Naukowa
677
Dodatkowymi elementami programu są dwa moduły graficzne. Pierw-
szy jest uruchamiany poprzez wybór z menu „Zanieczyszczenia” opcji „Dys-
persja” (F8) i przedstawia (rysunek 6) rozpraszanie smugi zanieczyszczeń, czyli
zmiany wartości współczynników dyspersji poziomej Sy [m] oraz pionowej Sz
[m] w funkcji zmiany odległości x [m] (równoległej do kierunku wiatru) od
receptora od emitera (przyjęto x
∈
<0 m; 10000 m>). Drugi moduł (rysunek 7)
uruchamiany przy pomocy polecenia „Stężenie” (F9) z menu „Zanieczyszcze-
nia” obrazuje graficznie rozkład stężenia Sx [µg/m
3
] substancji gazowej w osi
wiatru, na powierzchni terenu w odległości x [m] od emitora.
Rys. 7.
Rozkład stężenia zanieczyszczeń w osi wiatru i na powierzchni terenu
Rys. 7. Distribution of pollutant concentration in axis of wind and on the terrain surface
Uzupełniającym elementem programu jest możliwość generowania
raportu w postaci tabelarycznej, podzielonej na cztery kolumny – rysunek 8,
które zawierają informacje odnośnie nazwy parametru, symbolu, jednostki oraz
wartości parametru. Raport podzielony jest pod względem informacyjnym na
Jacek Piekarski, Marcin Lubierski
VI Ogólnopolska Konferencja Naukowa
678
dwie części. Część pierwsza zawiera dane początkowe (wejściowe), natomiast
w części drugiej przedstawione są wyniki obliczeń. Raport można wydrukować
lub zapisać w postaci pliku: tekstowego (*.txt), w formacie arkusza kalkulacyj-
nego Excel (*.csv) oraz w formacie html (*.html). Program umożliwia odczyt
oraz zapis wprowadzonych danych początkowych oraz współczynników obli-
czeniowych do pliku tekstowego (*.koc).
Rys. 8.
Raport symulacji komputerowej
Rys. 8. Report of computer simulation
Komputerowa symulacja rozprzestrzeniania zanieczyszcze
ń
w atmosferze
VI Ogólnopolska Konferencja Naukowa
679
3. Przykładowe obliczenia komputerowe
Dane wejściowe (Input data)
Wysokość emitera (od powierzchni terenu)
h = 44,0
[m]
Ś
rednica wewnętrzna wylotu przewodu emitujące-
go substancje
d = 1,6
[m]
Prędkość gazów na wylocie z emitera
v = 14,4
[m/s]
Temperatura gazów na wylocie z emitera
T = 493,80
[K]
Emisja NO
2
Eg = 43,9
[mg/s]
Prędkość wiatru (na wysokości anemometru)
Ua = 13,3
[m/s]
Ś
rednia temperatura powietrza dla okresu oblicze-
niowego
T0 = 270,13 [K]
Stan równowagi atmosfery
ATM = 4
[–]
Odległość receptora od emitera (równoległa do kie-
runku wiatru)
x = 500
[m]
Odległość receptora od emitera (prostopadła do
kierunku wiatru)
y = 500
[m]
Wysokość dla której określa się stężenie substancji
w terenie
z = 100
[m]
Współczynnik aerodynamicznej szorstkości terenu
Z0 = 0,18
[m]
Dane wyjściowe (Output data)
Stała obliczeniowa
m = 0,27
[–]
Stała obliczeniowa
a = 0,818
[–]
Stała obliczeniowa
b = 0,822
[–]
Prędkość wiatru na wysokości wylotu z emitera
Uh = 18,1
[m/s]
Emisja ciepła
Q = 4655
[kJ/s]
Zastosowana formuła
Hollanda
Wyniesienie gazów odlotowych
Dh = 2,6
[m]
Efektywna wysokość emitera
Hp = 46,6
[m]
Ś
rednia prędkość wiatru
(w warstwie od z = h do z = H
P
)
Usr = 18,3
[m/s]
Współczynnik obliczeniowy
A = 0,381
[–]
Współczynnik obliczeniowy
B = 0,218
[–]
Współczynnik dyspersji poziomej
(w warstwie od z = 0 do z = H
P
)
Sy = 61,47
[m]
Współczynnik dyspersji pionowej
(w warstwie od z = 0 do z = H
P
)
Sz = 36,06
[m]
Stężenie 1-godz. dwutlenku azotu
S = 2,5 E-16 [µg/m
3
]
Jacek Piekarski, Marcin Lubierski
VI Ogólnopolska Konferencja Naukowa
680
Wyniki przeprowadzonych obliczeń stężenia 1-godzinnego S [µg/m
3
]
dwutlenku azotu na powierzchni terenu dla przedstawionych danych wejścio-
wych i w odległościach x
∈
<0 m; 10000 m> oraz y
∈
<-1000 m; 1000 m>
(wartość ujemna składowej y oznacza odległość odmierzaną po lewej stronie
osi wiatru – do obliczeń podstawiana została wartość dodatnia) przy wyko-
rzystaniu programu „Atmo” zostały przedstawione w formie graficznej na
rysunku 9.
X [m]
Y [m]
0
-500
500
-1000
2000
10000
8000
6000
4000
0,01
0
0,05
0,04
0,03
0,02
S [ug/m
3
]
Rys. 9.
Stężenie NO
2
obliczone przy pomocy programu „Atmo”
Rys. 9. Concentration NO
2
calculated by programme „Atmo”
4. Podsumowanie
Istnieje wiele programów, przy pomocy których możemy symulować
rozprzestrzenianie się zanieczyszczeń w atmosferze. Programy te umożliwiają
modelowanie na dużych obszarach i z uwzględnieniem różnorodnych przemian,
którym podlega smuga zanieczyszczeń. Jednak są one kosztowne
Komputerowa symulacja rozprzestrzeniania zanieczyszcze
ń
w atmosferze
VI Ogólnopolska Konferencja Naukowa
681
i skomplikowane w obsłudze. Przedstawiony program komputerowy „Atmo”
jest aplikacją darmową (freeware), którą można pobrać ze strony internetowej
http://www.wbiis.tu.koszalin.pl/jacek. Służy do obliczania stężeń dla niewiel-
kich odległości (do 10 km). Pomaga wyeliminować czasochłonne obliczenia
przy wyznaczaniu stężenia zanieczyszczeń gazowych wg zaleceń Rozporządze-
nia Ministra Środowiska i jako taki może on stanowić pomocne narzędzie.
Literatura
1.
Lubierski M.: Numeryczne metody symulacji wybranych modeli rozprzestrzeniania
zanieczyszczeń w powietrzu atmosferycznym. Praca magisterska – promotor dr inż.
J. Piekarski, Politechnika Koszalińska, Koszalin 2002.
2.
Madany A., Bartochowska M.: Przegląd polskich modeli rozprzestrzeniania zanie-
czyszczeń atmosferycznych. Prace Naukowe PW, Inż. Środowiska 1995, Z. Nr 19
3.
Markiewicz M.T.: Przegląd modeli rozprzestrzeniania zanieczyszczeń w powietrzu
atmosferycznym. Prace Naukowe PW, Inż. Środowiska 1996, Z. Nr 21
4.
Nowicki M.: Parametry empiryczne w modelach dyfuzji zanieczyszczeń
w atmosferze. Zeszyt problemowy Nr X/84-85 Ochrona atmosfery z serii inżynier-
skie metody badania i obliczania stanu zanieczyszczenia powietrza atmosferyczne-
go, PZITS, Warszawa, 1984-85.
5.
Petterssen S.: Zarys meteorologii. PWN Warszawa 1966.
6.
Piecuch T. i inni: Spalanie i piroliza odpadów oraz ochrona powietrza przed szko-
dliwymi składnikami spalin. Wydawnictwo Politechniki Koszalińskiej, Koszalin
2002.
7.
Piecuch T. i inni: Ochrona powietrza. skrypt uczelniany Politechnika Częstochow-
ska, Częstochowa 1981.
8.
Rozporządzenie Ministra Środowiska z dnia 5 grudnia 2002 r. w sprawie wartości
odniesienia dla niektórych substancji w powietrzu, Dz. U. z 2003 r. Nr 1, poz. 12.
9.
Szklarczyk M.: Ochrona atmosfery. Wydawnictwo Uniwersytetu Warmińsko-
Mazurskiego, Olsztyn 2001.
10.
Woś A.: Meteorologia dla geografów. Wydawnictwo Naukowe PWN, Warszawa
1997.
11.
http://www.stronameteo.prv.pl/
12.
http://prognoza.icm.edu.pl/
13.
http://www.meteo.ids.pl/meu
14.
http://www.wunderground.com/global/PL.html
Jacek Piekarski, Marcin Lubierski
VI Ogólnopolska Konferencja Naukowa
682
Streszczenie
Zanieczyszczenia ulatujące ze źródeł emisji ulegają rozproszeniu
w atmosferze. Ich pomiar jest czasochłonny i kosztowny. Dlatego opracowano
matematyczne modele, które pozwalają z zadowalającą dokładnością określić
stężenia zanieczyszczeń. Jednym z takich modeli jest model gaussowski zale-
cany przez Rozporządzenie Ministra Środowiska z dnia 5 grudnia 2002 r. Auto-
rzy artykułu wykonali i przedstawili charakterystykę programu komputerowe-
go, opracowanego do obliczania stężenia zanieczyszczeń, wg wytycznych. Opi-
sali obsługę programu, zakres jego stosowalności oraz przykład. Przedstawiony
program jest prosty w obsłudze, zawiera szereg pomocniczych funkcji
i procedur ułatwiających pracę.
Computer Simulation of Pollutants Propagation
in the Atmosphere
Abstract
The pollutants outflowing from the emission sources undergo propaga-
tion in the atmosphere. Their measurement is time-consuming and expensive.
Therefore the mathematical models which allow to calculate concentration of
pollutants with satisfactory accuracy were worked out. One of such models is
Gaussian (recommended by Directive of Minister of the Environment). The
Authors of the paper created and introduced characteristics of computer pro-
gramme, worked out for calculation of pollutants concentrations, according to
guidelines. They described of programme maintenance, range of its applicabili-
ty as well as example. Introduced programme is easy in operation, contains
several auxiliary functions and procedures facilitating operation.