CC BY
24
Fayziev Khomitkhan, Doctor of Technical Sciences, Professor, Tashkent Institute of Architecture and Construction
Khozhiev T. K.,
Senior Lecturer, National University of Uzbekistan
Baymatov Sh. Kh., Senior Lecturer,
Tashkent Institute of Architecture and Construction
Rakhimov Sh.A., Senior Lecturer,
Tashkent Institute of Architecture and Construction E-mail: xomitxon@mail.ru
NUMERICAL SOLUTION OF THE BOUNDARY VALUE PROBLEM OF UNSTEADY FILTERING IN EARTH DAMS WITH ACCOUNT OF FILTRATION ANISOTROPY IN SOILS BY THE METHOD OF FINITE DIFFERENCES
Abstract. An algorithm for the solution of the problem of unsteady filtering in earth dams taking into account the anisotropy of soil using the finite difference method is given in the paper. The capabilities of the developed methodology, program for calculating unsteady filtering and the assessment of the reliability of the results obtained are illustrated by the example of filtration problems.
Keywords: earth dams, filtration, unsteady filtering, abrupt drop in water level in the upstream, rate of drop in water level, finite difference method, filtration anisotropy of soils.
During construction and operation of earth dams erected unsteady filtering in general case is described by three-dimen-
of clay and loamy sand soils, depending on the type of machinery used and production technology the soil of the dam body in some cases acquires anisotropic water permeability. Material of the dam body may be uniformly anisotropic or consist of alternating thin-layer isotropic layers with different values of filtration coefficient K. In the embankments erected by layered soil rolling or by means of hydro-mechanization, the value of Kx can exceed the value of Ky by the order of magnitude.
It is known that in vertical direction the permeability of loess soils (in natural occurrence) is greater than in horizontal direction. With other types of soils the opposite picture could be observed.
The study of filtration in earth dams with account of anisotropic water permeability has been studied by many researchers, including: N. N. Pavlovsky, V. I. Aravin, S. N. Numerov, P. Ya. Polubarinova-Kochina, V. P. Nedriga, L. N. Rasskazov, N. A. Aniskin, K. N. Anakhaev, R. A. Lyakhevich, M. Memari-anfard, Kh. Fayziev, S. N. Babakayev et al.
Solution of filtering problem in an anisotropic earth dam by exact methods is very difficult. Naturally, a qualitatively higher level of solving such problems is achieved using numerical methods of finite differences. Mathematical formulation of considered initial-boundary problem of anisotropic
sional nonlinear parabolic equations at various boundary conditions [2]. In particular, taking into account the symmetry along some spatial axes these equations can be reduced to two-dimensional equations in rectangular region Q = {a < x < b, c < y < d, tn < t <tk}, andwritten in the following form
Ar
dU dt
=4 ^
dx
A,
dU \ dx )
a d +A3—
dy
A,
dU dy
+ A,
dU dx
+ A,
dU dy
(1)
+ F,
at boundary conditions
dU
g 21
dx dU dx dU dy dU dy
x = a
x = b
y = c
+ g 12U(a,y,t) = q>i(y,t), c < y < d, (2)
+ g22U(b,y,t) = q>2(y,t), c < y < d, (3)
+ g 14U(x,c,t) = %(x,t), a < x < b, (4)
+ g24U (x ,d ,t ) = %(x ,t ), a < x < b, (5)
y = d and initial condition
U (x, y ,t )| t=n =y(x, y ),
(6)
13
NUMERICAL SOLUTION OF THE BOUNDARY VALUE PROBLEM OF UNSTEADY FILTERING IN EARTH DAMS WITH ACCOUNT..,
Where coefficients A,
(k = 0.6)
and function F can be
the functions of t, x,y,U, Ux and Uy, in particular, the constants used to describe equation (l); the coefficients gj are used to describe various boundary conditions. and W - are the given continuous functions. Assume that the given functions satisfy all the conditions under which there is a unique solution to the problem. To solve the problem (l) -(6), which in the general case is nonlinear, it is possible to apply various numerical methods of solution and linearization. According to the above, the main goal of the study is to create an algorithm and program for the numerical solution of this problem.
For numerical simulation ofthe problem (l) - (6) the method of variable directions (MVD) is used. First the continuous region Q is reduced to a discrete one QhT, with constant steps h = (b-a)/(m-1);I = (d-c)/(n-1);t =(tk -tn)/(p-1) along the Ox, Oy and Ot axes, respectively. After that, the approximation of the problem (l) - (6) at nodal points x{ = a + (i - 1)h, yj = c + (j -1)1 and ts = tn + (s -could be applied to the grid domain Q.hh and reduced to a finite-difference equation depending on the chosen method. Nonlinear function F on the right-hand side of equation (l) can be linearized using the Picard or Newton linearization method. Other methods of linearization also exist, but in practice the mentioned methods are the most
Aj
u 2 -U
11 ii
s+- s+- s+- s+-
A2ij h A2 ij
+(i -°)As6u
h
applicable. According to the Picard method, the values of nonlinear functions are determined from the values of the previous period of time and for the iteration case - from the previous iteration. In the algorithm, this is implemented by introducing a certain parameter p, so that its value can be used to choose a linearization method [3]. If the nonlinearity of F is related only to U, then the linearized function F has the following form:
F(t,x,y,U,UX,Uy)« F(ts,xt,yj,(U)j,(UX)j,(Uy)j)+ Fu(ts, x, yj ,(U )j ,(UX )j ,(Uy )j )(Uj+1 -Uj)
As a result of the linearization method used, at ¡ul = 0 the Picard method is obtained and at ¡ul = 1 - the Newton method. If the nonlinearity refers to Ux, Uy then they could be linearized in the same way: F (t, x, y ,U ,Ux ,Uy) *F(ts, x>, yj ,(U )j ,(Ux )j ,(Uy )j)+
+ n(FUx (ts, X,, yj ,mv ,(Ux )Sj ,(Uy )Sj )((Ux j-(Ux )Sj) + (8)
+ Fu (ts,xi,yj,(U)j,(Ux)j ,(Uy )j)((Uy j - (Uy )j))
For numerical solution of the problem it is advisable to apply absolutely stable implicit schemes, such as the MVD, with weight factor o, used to select calculation formulas using an explicit ( a = 0 ), implicit ( a = 1) and symmetric (a = 0.5) schemes. Bearing in mind the essence of the MVD, the calculation formulas are obtained first in x direction in the intermediate segment (s + l/2) in the form:
+ (l -
a)A3i, l [A
1 -us. us. -us..
',+1 - AS sJ-1
'4ij+1 j A4 ij
I
ii
s+- s
u 2-u 2
+ aAliiU+1j U-1j + 51j 2h
Uj+i -Uj-i 2I
and similar type of calculation formula in y variable in s + 1:
+ Fj +HiF'
i
us+2 -us. si si
(9)
A 2
A0 1,
us+1 -u
=(->A1if h
1
1
s+- s+
A21+1, Z
1
s+-
1 TT 2 TT 2
s+- u.2-u,.
- A 2 —ii-
A 1j h
s + M
+ aAs ijl
A 2
A4 1j+1
Us+11 -U1+1
- A 2
A4 1 j
1 us+1 -uj
+(1 j
1 1
1 +*< j
1 UÛ -us+1
2h
2I
11
ii 1 s+~ s ij-1 + f 2 + u F 2
ij r-1r Uij
1
U5+1 - U..2
(10)
Simplifying the coefficients standing before U^2 , Us
Uj+j in (9) and before Uj\, Uj1, U^ in (l0), the three-
diagonal systems of equations of special type are obtained
1 1 1
s+- s+- s+-
AjUi 2 + CjUn 2 + BjUi.2 = Dl, (ll)
solve the systems (ll) and (l2), the following formulas of the
scalar sweep method [2] are used:
B>s ds - AsBs _ _
, Pj = „ j . ( = 2,m-1, j = 2,n(l3)
a = —
Ci + Asah
C + AjaU
A.
ij i-1j u
ij i+1j
ij'
1
s+2 s+1
s+1
.. 1 + C+2U5.
1
+ Bj Ujh = Dj
Here
i = 2,m -1; j = 2,n -1;
as
. (12) Cj, bs , Dj
1111 A.2, Cj2, Bj2, Dj2 are determined from (9)-(l0). To
Uj = aU+p; (i = m - 1,m-2, ...,3,2; j = 2,n -1)(14) The values al and fif in (13) are determined by the boundary condition (2). In particular, for the boundary condition of the first kind obtained by g 11 = 0, g 12 = 1, (the Dirichlet problem) we get a^ = 0 and fi^ = Uj = (px. Con-
T
+
T
sidering this and other cases of boundary conditions on the left border x = a, the following algorithm is proposed. It is known that the three-point approximation of the boundary condition (2) of problems (l) - (6) has the form:
-U3 j 2 + 4U2 j 2 - 3U1j 2 1+1
3j 2j — + g 12U1 +2 =
2h ^. (l5)
Using this and formulas (ll) at i = 1, (14) at i = 2 the following system of equations is obtained relative to
1
1H— 2
11
— —
U 2 U 2 U
1j 2j 3j
1 1 1 A2jUj2 + c 2jU ;;2 + B2jU 3+2 =
Uj2 =«15U 2+2 +ß1
(16)
g 11
1
-U3j2 + 4U2j2 - 3U 1j2
+ 2hg 12Uj 2 = j. From the first equation we get
U;2 = D -Aipi;2 -c2jU2;2)/BSj, substituting it in the third equation (l6), the following could be determined
1
1+-Ui+2 =-
g 11 (c 2 j + 4B2 j )
g 11 ((j - 3B2j ) + 2hg 12B2j 2hB2^ 1j + gnD2j g 11 (A2j - 3B2j )+2hgBj '
-u 2+2 +
Hence:
= -
ßj =
g 11 (c 2j + 4B2j )
g 11 (A2j - 3B52j ) + 2hgB
2hB2 ^ 1j + g 11D2 j
(17)
1j ^ 11 (A2j - 3B2j)+2hgUB2'
In a similar form, to find the boundary value, the function
s+1 .
Um 2 is found from (3), on the right border x = b at 1 = m — 1 equations (ll), (14) and a three-point approximation in (3) are used:
1 1 1
A1 ,U1+2. + C1 ,U1+2. + B1 ,U1+2 = D
m-1 j m-2 j m-1 j m-1 j
11
1+- 1+—
U 2 =a1 U 2 + ß1
^ m-1j ^m-1^ m j H m-1
m-1j mj
1
m-1j >
(18)
1
1
1
3U 2 - 4U 1 + U 2. 1+-
mj m-1 j m-2 j T T -, 1
g21-j--j + g22Um j = 9 2j .
Determine Um -2 ■ from the first equation and substituting it in the last equation, using the second equation (l8) we get:
U 2
mj
2hKj12-g21K-1, +22,(4+ Cn )-1
g 21 (-1 - Bm-1 )+2hg 22Am-1 -gX -1 + Cm-1 )
(19)
The calculation method for solving the system of equations (12) in the direction of y variable is similar to the one in (13)-(14):
a
1
f-2
Bj2
1 11-
1+- 1+- 1+-
Cj2 + Va-2
11
2
ßf1 = D'j 2 A'j 2ßj-1 (i = 2,m -1, j = 2,n -1) (20)
- + - 1+- ^ '
1
1+- 1+— 1+— Cj2 + Vj2
1 1
Uj =ajS+2Uj11 + j2 (1 = 2,m-1; j = n-1,n-2,...,2). (2l)
Then the formula for the calculation of (20) at y = c can be determined:
/
1
1+-
a, 2 = —
11
1+- 1+-
C,22 + 4B 2
-! 2
V
/
ß1 2 =-
g 13 Ai 22 — 3g 13Bi 22 + 2g JB,
2lB* 2^3i + g 13Da 2
g 13 Ai 2
2 - 3g]3Bi22 + 2g 14/Bi2 and the function value at y = d :
1
1+2
(22)
(
4i"g23Din-l+g23
1 1A
1+- 1+-
4A,n-21 + C 2
n-1
u::
V
ßn—12
/
g23
1
3A»2 -1
+ WA»-' — g 23
1
-(23)
4Ain—2 + C»
The calculation process can be carried out by time or iterations. It is known that if the calculations are carried out by time intervals only, i.e. without iteration, the linearization is performed based on the previous time intervals. As for the iteration case, the calculations are carried out between the time intervals and such iteration process continues until the condition of convergence is satisfied
' "-<e (24)
max {|Uj
u j
(k+1)1 +1
-U(k
(k )1 +1
where £ is a sufficiently small number, 0 < £ < 1; k is the number of iteration, k = 0, l, 2
For this algorithm for solving problem (l) - (6) using the scalar sweep method, a program for calculating unsteady filtering in earth dams has been developed. Test problem. Consider the equation
d_H_ ~dt
= Kx -
x dx
H
dH dx
+ Ky dy
H
dH dy
(25)
initial conditions at K = Ky are taken the same as at isotropic filtering by the formula
NUMERICAL SOLUTION OF THE BOUNDARY VALUE PROBLEM OF UNSTEADY FILTERING IN EARTH DAMS WITH ACCOUNT...
H (t, x, y )|
\ >/ /\t=t =
0 = t H -
H1- H2
(26)
At anisotropic filtering, the initial position of the depression curve t = 0 is taken according to the formulas proposed in [4]. Boundary conditions are:
H(t,x,y)|x=a = H -fy (27)
H(t, x, y )| x= H2 (28)
dH
dy
dH dy
= 0
= 0
(29)
(30)
y=d
v,
Where H - is the pressure function in the computational domain, depending on coordinates and changing over time, m; t - is the time in days; a,b,c,d- are the borders of the domain of calculated profile; Kx, Ky - are the filtration coefficients in vertical and horizontal directions, respectively, they are constant, m/day; is the coefficient of water loss in soil of the dam body; H1,H2 - are the depths of the water in the upstream and downstream, respectively; v1 , v2 - are the rates of water drop in the upstream and downstream, respectively.
The developed methods allow solving the problems of unsteady filtration in earth dams without drainage and with various types of drainage systems (drainage prism, tubular drainage, layered drainage, etc.) with impermeable and permeable bases.
The possibilities of the developed methods of unsteady filtering and the assessment of the reliability of the results obtained are illustrated below using the filtration problems [3].
A homogeneous earth dam with a drainage prism in the downstream with an impermeable base is considered here. The initial depth of the upstream is 22 m; the laying of the upstream face is ml=3; the laying of the downstream face is m2 = 3. The ratios of soil filtration coefficients in horizontal and vertical directions is K/Ky = 4.0, the coefficient of water loss of soil of the dam body is y = 0.l; the rate of water level drawdown in the reservoir is 9 = 1.0 m/day, the level of the downstream is assumed unchanged and equal to H2 = 3 m.
Results of solving this non-stationary problem by the finite-difference method (FDM) and their comparison with the N. A. Aniskin and M. E. Memarianfard method [4] are presented in the table and (figure l).
Captions to the figure and table of the article "Numerical solution of the boundary value problem of unsteady filtering in earth dams with account of filtration anisotropy in soils by the method of finite differences".
aLB ' fuir [' \h
VHny MM
^ÏHE SM , 7 ►
-di
-id
-40
-10
It
3D
3)
50
Figurel. Results of calculations of unsteady filtering: 1 - initial depression curve at t = 0; 2 - positions of the depression surface by the N. A. Aniskin and M. E. Memarianfard method at t = 8 days; 3 - positions of the depression surface according to the FDM at t = 8 days
Table 1. - Results of calculations of unsteady filtering in a homogeneous dam with drainage prism and impermeable base
X, m FDM method FEM method Deviation%
-24 14 14 0
-20 13.96 13.75 1.5
-10 13.5 13.0 3.70
0 12.7 12.2 3.94
10 11.94 11.5 3.68
20 11.4 11.0 3.5
30 10.3 10.0 2.91
40 8.0 7.9 1.25
52 3.4 3.5 -2.94
In conclusion, it can be noted that the proposed numerical method for calculating unsteady groundwater filtration takes into account the anisotropy of soils in different directions. The calculations have shown that the re-
sults obtained using the FDM are quite true, qualitatively and quantitatively agree with the results of other methods compared, which indicates the reliability of numerical calculations.
References:
1. Aniskin N. A., Memarianfard M. E. Unsteady Anisotropic Filtering and its Effect on the Stability of Slopes of Earth Dam. Vestnik MGSU, 20ll.- No. 5.- P. 75-80.
2. Samarsky A. A. Theory of Difference Schemes.- M.: Nauka. l989.
3. Fayziev Kh., Khozhiev T. K., Khazhiev I. O. On One Algorithm and Computational Experiment for One Problem of Earth Dam. Proceedings of the International Conference "Actual problems of applied mathematics and information technology -Al-Khorezmi 20l6". - November 9-l0. 20l6.- Bukhara.- P. 47-49.
4. Anakhaev K. N., Lyakhevich R. A. To the Filtration Calculation of Anisotropic Earth Dams with layered drainage // Hydro-technical Construction, 2006.- No. 9.- P. l9-22.
