CC BY
21
УДК 519.6
M. M. BILIAIEV, P. S. KIRICHENKO (DIIT, Dnipropetrovs'k, Ukraine)
NUMERICAL SIMULATION OF THE BLACK SEA POLLUTION FOR THE CASE OF THE MINE WATERS DISCHARGE
На 6a3i розроблено! 3D чисельно! моделi виконано розрахунок забруднення aKBaTopii' моря. Модель базу-еться на iнтегруваннi рiвняння конвективно-дифузшного переносу домiшки та моделi потенцшного руху. Наводяться результати обчислювального експерименту.
Ключовi слова: забруднення акватори моря, шахтовий водоскид, рiвняння конвективно-дифузiйного переносу домшки, модель потенцiйного руху, тривимiрна чисельна модель
На базе построенной 3D численной модели выполнен расчет загрязнения моря. Модель основывается на численном интегрировании уравнения конвективно-диффузионного переноса примеси и модели потенциального течения . Приводятся результаты вычислительного эксперимента.
Ключевые слова: загрязнение акватории моря, сброс шахтных вод, уравнение конвективно-диффузионного переноса примеси, модель потенциального движения, трехмерная численная модель
The 3D numerical model was developed and used to simulate the sea pollution after mine waters discharge. The model is based on the numerical integration of the K-gradient transport model and the model of potential flow. The results of numerical experiment are presented.
Keywords: sea pollution, mine waters discharge, K-gradient transport model, model of potential flow, 3D numerical model
Introduction
The mining the ore is the leading branch of Ukraine economy. But this activity results in a huge amount of mine waters. Utilization of mine waters in Ukraine is a problem of great importance. The mineralization of mine waters covers the range 19... 36 g/l. The annual total amount of mine waters in Kryvyi Rih region (the central part of Ukraine) is about 21-106 m3. Now these mine waters are collected in special ponds (Fig. 1) and then discharged in rivers.
Fig. 1. Mine waters discharge into the pond Svistunova
The mining the ore is planned to be increased that is why the problem of utilization of mine waters attracts again the attention. In 1988 the project of mine waters discharge in the Black Sea was de-
veloped in the USSR but the ecological study was not carried out at that time. Now this project is under consideration and the problem of the ecological study is appeared again. The engineers working on the correction of the old project to meet the modern demands need the tool to predict the sea pollution for the case of mine waters discharge into the Black Sea. The Regulation Instructions based on the analytical model and which are used now in Ukraine to predict the sea pollution could not provide the solution of the problem considered This paper presents a numerical model to simulate the flow field and pollutant dispersion in vicinity of the discharge pipe opening which is situated in the sea. The hydrodynamic model of inviscid flow is used to predict the velocity flow field which is formed as a result of the sea flow and mine waters jet interaction. This model does not consume much of computing time.
Mathematical model of pollutant dispersion
To simulate the pollutant transport in the sea the gradient transport model is used [1, 2, 3]
dC
dt +
__d_ dx
duC dvC д(w - wg )C
dx dy
V x
dC
dx ) dy \ dy
V y
dz
dC
cC _
dz ^ z dz ^
h^ Q {t)ô(x - xt )S(y - yt )5(z - zt ), (1)
© Biliaiev M. M., Kirichenko P. S., 2011
where C is concentration; u, v, w are the velocity components; c is the parameter taking into account the process of pollutant decay; ^ = , ^ , ^z) are the diffusion coefficients; Q is the intensity of pollutant emission.
To simulate the mine waters gravity fall as the result of the density difference between them and sea water the parameter wg is used. This parameter is determined in experimental studies.
The transport equation is used with the following boundary conditions:
• at the inlet boundary: Ciniet = CE, where CE is the known concentration;
• at the outlet boundary: in numerical model the condition C (i +1, j, k) = C (i, j, k) is
used (this boundary condition means that we neglect the process of diffusion on this plane);
• at the top boundary and the bottom
surface: dC = 0 .
dn
In the numerical model the following approximations for sea averaged speed is used:
um = a0 + a1 V + a2 H + a3 V2 + a4 H2 + a5 VH + + a6 V2H + a7 VH2,
where ai (i = 1, 2, ..., 7) are the coefficients, which are shown in Table 1; H is the sea depth; V is the wind speed.
This approximation is used if the following conditions are fulfilled:
2.0 < V < 20(m/s), 1.5 < H < 50(m).
Table 1
Coefficients to calculate the averaged sea speed and parameters of diffusion
V < 6 m/s
ûi Ci
a0 = 3,613-10-2 c0 = 599-10-4
a1 = -2,751-10-3 Cj = 5,347-10
a2 = 1,108-10-2 c2 = -3,681-10-4
a3 = 1,461-10-3 c3 = -1,469-10-4
a4 = 9,729-10-6 c4 = 5,669-10-6
a5 = -7,189-103 c5 = 1,426-10-4
a6 = 9,925-10-4 c6 = 2,276-10-6
a7 = -3,875-10-6 c7 = -2,401-10-6
The diffusion coefficients are calculated using the following approximations:
^x = 0.032 + 21.8u2m ;
= c0 + c1V + c2 H + c3V2 + c4 H2 + c5 VH +
+ c6 V2 H + c7 VH2,
where ct (i = 1, 2, 7) are the coefficients shown in Table 1.
It is worth noting that these approximations are recommended in the State Regulations, which are used in Ukraine to predict the pollution of the sea if the problem of waste waters injection into the sea is under consideration. The State Regulations do not indicate how the diffusion coefficient along Y axis must be calculated so in the model developed we consider that ^y = .
Numerical model
The calculation of pollutant dispersion is carried out on the rectangular grid. Main features of the finite difference scheme for transport equation are presented below.
The time dependent derivative is approximated as follows:
d C d t
n + 1 n
C ijk - C ijk
A t
At the first step convective derivatives are represented in the following way:
duC du+C du - C
dx ---1 dx i--. dx
cvC dv+C dv'C
cy cy dy
dwC dw+C dw~ C i--
where u + =-
dz
u + lui
dz
dz
u — u
v+v
2
■; u =-
v =-
v =
v —Iv + w + w
w =•
2
w — |w|
; w =■
2 2 2 At the second step the convective derivatives are approximated as follows:
du+C ^ u++1, j, k Cj — uijk C"—+\, j ,k dx Ax
du C ^ ui+1, j,k Ci++i) j,k — uijk C"jk d x Ax
dv+C v+i,j+1,kCijk — v++kCi, j—1,k
= L+ C
= L— Cn
n+1
cy
Ay
= L+Cn
2
dv C vi, j+1,kCi, j+1,k VijkCijk _ L- cn+1
_ LyC '
dy
Ay
cw+C ^ j k+iCl}k _ wl]kCl, jk-i _ L+cn+i.
dz
Az
"'z '
as
_ / _ X /-»n+1 /-»n+1 /-»n+1 /-»n+1
_( _ \ ~ ~ Ci+1'j• k _ Cijk _ ~ Ci-j• k _ Ci'-1.-
dx v x dx J x 2 x
A x
A x2
_ m:C"+1 + MiCn+1
f nn+1 nn+1 nn+1 nn+1
d \ dC \ ~ Ci, j+1,k _ Cijk _ Ci, j,k _ Ci, j_1,k
-»n+1 . r+ n+1
_ M- Cn+1 + M+ C ;
d f dC ^ _ Cn +1+1 _ C,
n+1 _ n+1 n+1 _n+1
i, j,k+1 _ Cijk ~ Ci, j,k~ Cij,k-1
: (M+ Cn+1 + M- Cn+1 + M+ Cn+1 +
\ xx xx yy
+ ^ Cj _ 4 ( M+x Ck + M-JCn + M+y Ck + + My Cn + M+ Ck + Mz Cn ) ; (3)
at the second step ( k _ n + -1; c _ n +1 )
2 4
the difference equation is
ik /nc 1
ijk _ 'Jk +1 (l Ck + L_ Ck + LL Ck) +
2 \ x y z )
At
dw C " k+1Ci,j,k+1 _wykCi,j,k _ l_C"+1 dz Az z
The second order derivatives are approximated
-»n+1 /—»n+1 /—»n+1 /-»n+1
+ - ck =1 (mxck + M+Cc + M- ck +
4 iJ 4 \ x x yy
+ M+y Cc + M- Ck + M+ Cc) ; (4)
yy zz zz f ' v '
3 1
• at the third step (k = n + —; c = n + —) the
4 2
expression (4) is used;
3
• at the fourth step (k = n +1; c = n + —) the
expression (3) is used. At the fifth step (if we want to take into account the pollution source influence) the following approximation is used:
fel^f ^z ~Az2 ' Az2
_ M-Cn+1 + M+7Cn+1 .
zz zz
In these expressions L+, L_, L+, L_, L+, L_,
M+x, M-Ç, ... are the difference operators. Using these expressions the difference scheme for the transport equation can be written as follows:
Cn+1 Cn
iJk _ iJk + L+ Cn+1 + L Cn+1 + L+ Cn+1 +
At x x y
T- /"»n+1 . T+ n+1 . T _ S~in+1 . _/"»n+1
-LyC + Lz C + Lz C +a C1Jk _
5 n+1 5 n
Ci, j ,k_ Ci, j,k _ N q i(t"+12)
At A x A y A z
^ .
Function Si is equal to zero in all cells accept the ones where the «1» source of pollution is situated.
This difference scheme is implicit and absolutely steady but the unknown concentration C is calculated using the explicit formulae at each step (so called «method of running calculation»).
To simulate the flow field the model of potential flow is used. In this case the governing equation is
d2P d2 P d2 P
dx2 dy2 dz2
+ M- Cn+1 + M+ Cn+1 + M'Cn+1) .
yy zz zz f
Solution of the transport equation in finite-difference form is split into four steps during the time step of integration dt:
• at the first step (k = 4) the difference equation is:
C n+k c n i
l] ~ l] +1 (l +xck + L\ C k + L+z C k ) + At 2V x y '
_ 0.
(5).
where P is the potential of velocity.
The components of velocity are calculated as follows
dP dP dP u _—; v _— ; w _ —. dx dy dz
The boundary conditions are as follows:
dP
- at the bottom and the upper boundary: — _ 0,
dn
where n is normal unit vector to the boundary;
dP
- at the entrance boundary: — _ Vn, where Vn is
dn
the known speed;
- at the exit plane:
P _ P * ( x _ const, y) + const .
To solve the Laplace equation for potential of velocity the Libman's method is used. In this case the approximation of this equation is as follows:
P - 2P + P P - 2P + P
Ti+1, j,k j,k i-1, j,kri, j+1,k i, j,k i, j-1,k
Ax2
P -2P +P
ri, j ,k+1 j ,k i, j ,k-1
AZ2
Ay2
= 0.
From this expression we can find the unknown function Pij
=a
i (P - P P - P
1 1 i+1, j, k 1 i-1, j, k 1 i, j+1, k 1 i, j-1, k
Ax2
Ay2
+
P - P
1 i, j,k+1 1 i, j,k-1
AZ2
r
where A =
Ax2 Ay2 Az2
The components of velocity vector are calculated on the sides of computational cell as follows:
4 ■ ■ i
i,} ,k
P - P
i, j ,k
Ax
Pi, j,k - Pi, j-1,k
Ay
Pi, j ,k -P ri, j ,k-1
Az
A code based on this numerical model was developed. This code was used to predict the sea pollution after mine waters discharge into the Black Sea near settlement Zaliznyi Port («Iron Port») that is situated in Kherson region of Ukraine.
Results
A numerical experiment was carried out for the following initial data. The wind speed equals 3 m/s; the mine waters speed at the pipe opening is 1.27 m/s; the sea depth is 12 m; the pipe diameter is 1 m; the concentration of pollutant in mine waters is 100 units (dimensionless value). The dimensions of the computational region are 12 x 16 x x 12 m. The numerical experiment was carried out for two values of wg: wg = 0.01 m/s, wg = 0.03 m/s.
The results of the numerical experiment are shown in Figs. 2-5.
It's clear that the increase of wg value causes the decrease of polluted area. So the experiments to estimate this value must be carried out with the needed accuracy. It's obvious that the most intensive polluted area is formed near the opening and
is about 11 m long and at this distance the concentration decreases intensively.
—8.1031'133 u 0.967E+B2 r —B.'JlM'Lt: J
0.833E+0Z i _0.7B7E«0Z n 0.700E+0Z a t
0.567E+0Z e J.1VII.FE 0.433E+02 —0.3G7E+02 y
0.300E+02 —0.233E+02 0.167E+02 —0.1B0E+02 0.333E+01
coordinate x
Fig. 2. Pollutant concentration near the pipe opening; wg = 0.01 m/s (side view, section y = 5.5 m)
Fig. 3. Pollutant concentration near the pipe opening; wg = 0.03 m/s (side view, section y = 5.5 m)
Fig. 4. Pollutant concentration near the pipe opening; wg = 0.01 m/s (top view, section z = 1.5 m)
Fig. 5. Pollutant concentration near the pipe opening; wg = 0.03 m/s (top view, section z = 1.5 m)
Conclusion
The Ukrainian Environment Protection Law demands that the concentration of pollutant must not exceed the permitted concentration level at the distance 500 m from the point of discharge. From this point of view the discharge considered meets this demand.
REFERENCES
1. Израэль, Ю. А. Экология и контроль состояния природной среды [Текст] / Ю. А. Израэль. - М.: Гидрометеоиздат, 1984. - 560 с. (in Russian: Izrael', Yu. A. Ecology and Monitoring the Natural Environment Condition [Text]. - Moscow: Gidro-meteoizdat, 1984. - 560 pp.)
2. Озмидов, Р. В. Горизонтальная турбулентность и турбулентный обмен в океане [Текст] / Р. В. Озмидов. - М.: Наука, 1968. - 204 с. (in
Russian: Ozmidov, R. V. Horizontal Turbulence and Turbulent Exchange in the Ocean [Text]. -Moscow: Nauka, 1968. - 204 pp.)
3. Саркисян, А. С. Численный анализ и прогноз морских течений [Текст] / А. С. Саркисян. - Л.: Гидрометеоиздат, 1977. - 182 с. (in Russian: Sarkisyan, A. S. Numerical Analysis and Forecast of Sea Flows [Text]. - Leningrad: Gidrometeoiz-dat, 1977. - 182pp.)
Received by Editorial Board: Nov. 11, 2010.
Accepted for publication: Nov. 23, 2010.
