УДК 539.3
Modeling of the Dynamics of Spreading of Spilled Hydrocarbons Taking into Account the Gravity-Capillary Interaction
Natalia V. Molokova* Margarita A. Konnykh
Institute of Space and Information Technology, Siberian Federal University, Kirensky, 26, Krasnoyarsk, 660074,
Russia
Received 10.05.2012, received in revised form 10.06.2012, accepted 20.06.2012 It is covered a two-dimensional model of filtering, taking into account the movement of air and hydrocarbon pollutants into the porous soil. The model includes a system of partial differential equations with additional conditions. Among the differential equations are mass balance equation element in a porous medium-inseparability equation, and the differential equations of motion. For circuit system entered the equation of state before pollutant and the environment. Initial and boundary conditions correspond to the filtration process, beginning with the ground surface and initial stage of oil products ladling. A comparative analysis of the results of mathematical modeling with experiments.
Keywords: two-dimenstonal model of filtering, geofiltration, mathematical modeling, porous medium, hydrocarbon pollutants.
The relevance of this work is determined by the need to promptly receive and process information about intensity and nature of the spread of spilled hydrocarbons in order to develop an optimal system of measures for the prevention and elimination of oil contamination. To increase significantly the efficiency of process of developing an environmental action is possible with the help of computer modeling and information technologies.
The analysis of the literature shows that at the present time there is a wide range of models of filtration of hydrocarbons in porous environments and software systems that implement these models. However, the use of the proposed models for the majority of contaminated soil is limited for several reasons: first, difficulties connected with the subsequent model equipment by the adequate initial data; secondly, the slow account of tasks; and thirdly, the high cost of the program complexes.
Therefore it was decided to create a new information system, which would eliminate these shortcomings [1]. In this connection object of research is process of petropollution of soil as a result of emergency leaks and spills, and an object of research is applied aspects of modeling of a filtration of liquid hydrocarbons in the porous environment.
The purpose of the work was to model the dynamics of the spread of spilled hydrocarbons, taking into account relevant factors to identify areas with a high degree of petropollution, where carrying out of regenerative and cleaning works is necessary.
In order to achieve the objective it was necessary to solve the following tasks:
* [email protected] © Siberian Federal University. All rights reserved
1. To analyse a research problem in the scientific domestic and foreign literature.
2. To construct mathematical model for an estimation of intensity and character of the spatial-temporal distribution of petropollution in soil environments. Thus the model should be adequate to natural conditions.
3. To develop the program complex, which allows to conduct research operatively and to receive trustworthy information about intensity of vertical and horizontal migration of the contaminant in soil, about degree of petropollution and the area of its distribution.
4. To solve a Geofiltration problem and to carry out the comparative analysis of the received results with experimental data.
1. The Analysis of Experimental Data
As in any method of scientific knowledge, in modeling there are experimental and theoretical parts.
The initial information for a substantiation of system of the basic equations of model has been received from the analysis of a series of experiments: the analysis of the mechanism of pollution well permeable and less permeable environments; determination of hydrocarbon migration dynamics in macrohomogeneous air-dry soil; analysis of the mechanism of spreading of a pollutant on the boundary between two environments of different permeability [2].
Analysis of experimental data has shown, that getting on the surface of the earth, liquid hydrocarbons are beginning to seep in pores of aeration zone, dominated by vertical migration. When petroleum products meet on their way less permeable layer there is their accumulation and spreading in the horizontal direction. Migration processes of hydrocarbons in soils are defined by their properties, parameters of environment and are controlled by laws of the theory of a filtration.
2. Conceptual Statement of a Problem
1. The incompressible liquid (p = const, p = const) in not deformable ground (m = const, k = const).
2. We believe that the petropollutant spreading occurs along an axis z, directed from a surface the considered area vertically downwards. Lateral borders of area don't influence filtration process.
3. In the first approximation soil saturation only by air is considered, pressure of air is considered equal to the atmospheric.
ni
4. The sum of saturations is equal to unit si = 1 then for a diphasic filtration si + s2 = 1,
i= i
therefore from two saturations only one is independent and the notation is introduced s = s1 — oil saturation, s2 = 1 — s.
5. At a diphasic flow of immiscible liquid, pressures in each of phases aren't equal among themselves (p1 = p2), pc(s) = p2 — p1, indexes 1, 2 concern a pollutant and gas accordingly.
According to the accepted assumptions the set of parameters has been defined.
Input modeling parameters are: properties of a pollutant (density pi, viscosity pi); properties of the second phase (density p2, viscosity p2); characteristics of the soil (porosity m; permeability k); relative phase permeability (k1(s),k2(s)); depth and width of filtration area (L1, L2); limiting values of a saturation s* = 0, 9; s* =0,1; a saturation on the top border of area S0 = 1,0.
Parameters of influence of environment are: acceleration of free falling g; capillary pressure
Pc(s).
Output parameters of the model: the degree of petropollution, which is determined by the oil saturation s1; the depth and width of its dissemination; the speed of formation of a zone of petropollution.
3. Mathematical Model
On the basis of the laws controlling a current of hydrocarbons, the full balance equations of a diphasic filtration taking into account operating factors (interaction of capillary and gravitational forces) have been written out.
m— = div dt
, ki(s)
k-(VPi - Pig)
Pi
P2 = Pi + pc(s), si + S2 = 1, for s = si value S2 = 1 — s,
0(1 — s)
dt
div
kk2(s)(vp2 — p2g) P2
(1)
(2)
(3)
(4)
From the equations of continuity and Darcy's law the equation for determining the pressure as a function of coordinates and time was derived.
div
k ( ki(s) + k2(s)
Pi
P2
V Pi
—div
kMf! v pc(s) P2
dki(s) „ dk2(s) _ + k-—Pig + k-— P2g.
pidz
P2$z
(5)
The operator of the model includes system of the differential nonlinear equations of parabolic and elliptic types. This model describes the process of unsteady filtration of hydrocarbons in different soil environments with variable filtration coefficients.
Filtration coefficient includes the important parameter — the relative phase permeabilities k1 (s) and k2(s), which, like the pc(s), are experimentally measuring functions of saturation, forcing out phase. Their typical form and the empirical formulas for a system of oil — gas are shown in Fig. 1.
ki(s)
(^)ai ,s* < s < 1 s*
0, 0 < s < s*
(6)
k2 (s)
-)°2, 0 << s < s*
0, s* < s < 1
(7)
Let us formulate the problem, that is, we define the conditions at the initial time and the boundary conditions at the boundaries of the study area. The system of equations (1-7)
*
(
*
Fig. 1. The relative phase permeabilities
Fig. 2. Formal representation of subject domain (Qt — area of a filtration of a pollutant, Q0 — the area of oil spills, Q — computational area)
is considered in the area Q x (0; T], where — Q is three-dimensional cube with boundary r, Q = [0; li] x [0; I2] x [0; I3] — are shown in fig. 2.
The initial-boundary conditions: at t = 0: s = 0, (x, y, z) E Q, p = Patm, (x, y, z) E Q, Patm — atmospheric pressure. The boundary conditions were chosen in accordance with the features of the considered surface.
The border with the atmosphere (r1). On the surface of contact with the atmosphere (r), s = 1, (x,y,z) E Qo, p = pigz + Patm, (x, y, z) E Qo, s = 0, (x,y,z) E Q \ Qo, P = Patm, (x,y,z) E Q \ Q0. On the border (r2) between the phases the oil pressure is given by P = Patm + Pc (s).
Impermeable boundary (Q3). Fluid cannot pass through an impenetrable boundary, so the
k1 (s)
condition of impermeability can be written as k-(VP1 — p1g) • n|r =0, n — normal vector
13
to the boundary of the r3.
A complex mathematical model, describing the intensity and the nature of the spatial-temporal distribution of oil pollutants in the upper part of the aeration zone is subject to further
the numerical solution.
4. Approximation
In the numerical solution of the reduced system of equations was used finite-difference approximation. To construct a difference scheme for the original problem, we choose a uniform grid with step on the variable z, with step on the variable x, with step hy on the variable y and step t on the variable of the time t. Get the grid area:
Wh,T = Zii(0 < ¿1 ),Xi2 = (0 < ¿2 < n®),yi3 = ¿3hy(0 < ¿3 < ny), tj = jT, 0 < j < nt, where nz = L1/hz,nx = L2/hx,ny = L3/hy, nt = t/t. We replace the derivatives at the interior nodes Wh,T by finite-difference relations taking into account the representation of the variable coefficient and we will get the difference scheme:
mst = (^Aa(S)Pg „ - G ^ (8)
(sj+11 + sj+1A
where Aa±Q,5 (s) = k 8 a±12 , at s = 1 phase permeability k2(s) addresses in a zero,
then the difference approximation of the equation for the pressure we can write in the following form:
Pxa x a0, (9)
at 0 < S < 1, PX+1 = P0tm + Pcj+1(sj+1) (10)
s0,»,», = 0, Piii2i3 = Patm, j = 0, ¿1 = 0, ...,nz - 1; ¿2 = 0, ...,nI - 1; ¿3 = 0, ...,ny - 1;
il«2»3
cj - I
SQi,i, = 1, ¿1 = 0; 0 < ¿2, ¿3 < 3; sQi2i, = 0, ¿1 = 0; 3 < ¿2, ¿3 < nKy - 1; pQi2i3 = P1g + Patm, ¿1 = 0; 0 < ¿2 ¿3 < 3; pQ»2. = Patm, ¿1 = 0; 3 < ¿2, ¿3 << nKy - 1; (11)
sj+1 - sj 2 • SN sN = - 2 AJ
N-1/2pN - P1g j = 1, nt-1
at 0 < s < 1, PJ+1 = Patm + Pj+1(sj+1) (12)
sQ1i2i3 = 0, j = 0, ¿1 = 0,..., n2 - 1; ¿2 = 0, ...,n, - 1; ¿3 = 0, ...,ny - 1; sQi2i, = 1, ¿1 = 0; 0 < ¿2, ¿3 < 3; sQi2i, = 0, ¿1 = 0; 3 < ¿2, ¿3 < nKy - 1;
pQi2i3 = P1g + Patm, ¿1 = 0; 0 < ¿2 ¿3 < 3; pQi2i, = Patm, ¿1 = 0; 3 < ¿2, ¿3 << nxy - 1; sj+1 - sj 2 .
sN sN 2 \j j -1 -i
—T— = -hAN-1/2PN - P1g j =1nt-1.
The formulated statement of a problem where the unknown functions are P, S allows to result the scheme of the numerical decision [3].
1. We enter the initial data (including information from experiments).
2. According to the known on the j-th time step Pj, and the boundary conditions for the P we define P from the difference ratio (9) of the analogue of the equation (5) with using the iterative Gauss-Seidel method or (10).
3. According to the found Pj+1, known Sj, boundary conditions for Sj, define Sj+1 from grid (8) analogue of equation. Equation (1) has the parabolic type. It is expedient to apply the iterative implicit finite difference schemes.
4. Then, the value P|+1 = P|+1(Sj+1) is determined by the found Sj+1. Account values on the (n + 1)-th step is completed, proceed to the next step and etc.
5. Results of calculation of the unknown functions P, S are displayed in tabular and graphical forms.
5. The Complex of Programs
To implement the proposed computational algorithm was developed software complex. As a tool for its creation was chosen technology of visual programming Borland Delphi, which has allowed to eliminate effectively the disadvantages of other complexes. The program complex includes three modules and a dialogue window. Program work begins with loading of the initial data in editing fields. Step-by-step calculation of a task is made by control program of the complex. The control program interacts with calculation modules through the core of the complex, which carries out work with the data using a set of appropriate interfaces. The first module is intended for an estimation of intensity of pollution of soil environments on depth. The second module is intended for an estimation of speed of formation of a zone of pollution and its sizes. The multimedia module is intended for an estimation of area and degree of pollution of soil environments at a filtration of various hydrocarbonic pollutants in a dynamic mode.
6. Results of Numerical Calculations
According to the proposed scheme series of calculations has been carried out [3].
Figure 3 shows the numerical calculations of the dynamics of saturation over the depth and width, which show ambiguous character of spatial-temporal distribution of oil pollution. For example, if environment is presented by loamy soils the horizontal filtration prevails over the vertical. The petropollution leaks down on depth of 0,7 m, and spreads on a horizontal plane on distance to 1 m for 37 hours. According to experimental data spilled hydrocarbons are fixed in the upper layers of the loamy and clay soils. Qualitatively it is explained by low bandwidth of the environment. And on the contrary, if environment is presented by less disperse breeds (in this case sand) vertical migration prevails. The petropollutant extends deep into to 1 m, and across to 0,6 m for 24 hours. It was also identified three zones with various degree of pollution: strong (value of a saturation is allocated by black color), medium (value of a saturation is allocated by dark gray color) and weak (value of a saturation is allocated by light gray color).
On fig. 4 the result of calculation of dynamics of distribution of a pollutant on border of section of two environments with various permeability is presented.
Contaminant is distributed almost in the vertical direction. Above an impermeable boundary of two environments, which is on depth of 1 m, there is its accumulation and horizontal distribution. Explanation to that is, that the gravitational pressure appears more strongly capillary counteraction during all estimated time.
In general, the process is proceeding quite slowly. But it should be noted that the maximum speed of filtering were identified for species with a minimum content of the sand fraction, and the minimum speed for the rocks with a maximum content of sandy fraction. Filtration time
in the sand was 22 h, the speed of advancement of the front is approximately equal to 5 • 10-5; filtration time in a clay sand accounted for 28 h, the speed of progress of the front approximately equal to 2, 5 • 10-5; filtering in sandy loam amounted to 33 hours, the speed of advancement of the front is approximately equal to 1, 67 • 10-5; and the time of filtration in the loam to 39 hours, the speed of advancement of the front is approximately equal to 1, 25 • 10-5, which is consistent with experimental data. Consequently, a significant amount of petropollutant migrates in the direction of improving the filtration characteristics of the soil, and here we should expect the maximum pollution.
Calculations also show that the gravity-capillary interaction significantly affects the dynamics
Fig. 3. Dynamics of a petrosaturation on depth z and width x in environments: a) loam (m = 0, 25; k = 10 - 15 m2; T = 37 h); b) sand (m = 0, 25; k = 10-12 m2; T = 24 h)
Fig. 4. Petropollution distribution to a vertical cut during the various moments of time T(daity) —• a) 1,48; b) 1,85
of the current depth of a petropollutant in soil. The increasing of the gravitational parameter leads to a reduction of s*, and then the decisions take the form of a similar moving up stairs. It tells about increase in a front petrosaturation, i.e. about full penetration of a pollutant deep into soils. On the other hand, if more strongly capillary counteraction the orientation effect
decreases — the front form becomes more smoothed. An explanation to that is, that capillary forces impede the movement of a pollutant into the depth. This fact should be considered at calculation of a zone of pollution and a choice of ways of clearing.
7. Check of Model Adequacy
Results of comparison of calculated values with experimental data allow to draw a conclusion that authenticity of model has reached 82% since the size of the minimum root-mean-square deviation has made 18% at the chosen scale factors and the main parameters of the process (duration, depth and width of distribution of a pollutant). Thus, the offered model is adequate to real process and can be used for an estimation of degree of oil pollution and its sizes as a result of emergency leaks and spills.
Conclusion
1. The three-dimensional mathematical model for research of dynamics of distribution of the spread of spilled hydrocarbons, taking into account the gravitational-capillary interactions was offered.
2. The effective computing algorithm that takes into account the specific features of a particular mathematical tasksis constructed. The special feature of this algorithm is that the calculations are carried out only in the field of distribution of hydrocarbons, which allows to considerably reduce the computational cost with the same accuracy.
3. The specialized program complex which possesses a number of advantages is developed: allows to conduct researches operatively; gives the chance to solve problems for which there are no standard computing algorithms.
Further at calculations it is planned to take into account such physical characteristics of objects of research as humidity of soils and ambient temperature.
References
[1] N.V.Molokova, Mathematical modeling of pollution distribution in the porous environment, Supplement to the journal "Open education". Materials of the All-Russia conference "IV Siberian congress of women-mathematicians", Krasnoyarsk, (2006), 130-138 (Russian).
[2] N.V.Molokova, The decision of geofiltrational problems by means of mathematical modeling, Vestnik of the Siberian State Aerospace University, (2008), no. 3(20), Krasnoyarsk, Siberian State Aerospace University, 31-34 (Russian).
[3] N.V.Molokova, Numerical research and the analysis of process of pollution of the porous environment, Vestnik of the Siberian State Aerospace University, (2008), no. 4 (21), Krasnoyarsk, 20-23 (Russian).
[4] N.V.Molokova, Mathematical modeling of petropollution processes in the porous environment, Vestnik of the Siberian State Aerospace University, (2010), no. 5 (21), Krasnoyarsk, 142-148 (Russian).
Моделирование динамики распространения разлитых углеводородов с учетом гравитационно-капиллярного взаимодействия
Наталья В. Молокова Маргарита А. Конных
Рассматривается математическая модель двухфазной фильтрации, учитывающая движение углеводородных загрязнителей и воздуха в пористом грунте. Модель включает в себя систему уравнений в частных производных с дополнительными условиями. В число дифференциальных уравнений входит уравнение баланса массы в элементе пористой среды - уравнение неразрывности, а также дифференциальные уравнения движения. Для замыкания системы вводятся уравнения состояния рассматриваемого загрязнителя и среды. Начальные и граничные условия соответствуют фильтрационному процессу, начиная с поверхности грунта и начальной стадии разлива загрязнителя. Проводится сравнительный анализ результатов математического моделирования с экспериментами.
Ключевые слова: модель двухфазной фильтрации, геофильтрация, нефтезагрязнение, пористая среда, математическое моделирование.