CC BY
37
УДК 532. 546
ИДЕНТИФИКАЦИЯ КОЭФФИЦИЕНТА ФИЛЬТРАЦИИ ПОЧВОГРУНТОВ МЕТОДОМ РЕГУЛЯРИЗАЦИИ
М.У. Мурзакматов, Б.А. Байболотов
Иссык-Кульский государственный университет им. К. Тыныстанова, Республика Кыргызстан
Представлена членом редколлегии профессором Г.М. Куликовым
Ключевые слова и фразы: коэффициент фильтрации; метод конечных элементов; уровень грунтовых вод.
Аннотация: Рассматривается применение метода регуляризации А.Н. Тихонова к идентификации коэффициента фильтрации в двумерном уравнении Бус-синеска.
Коэффициент фильтрации является одним из основных гидрогеологических параметров, характеризующих водоносные пласты. Для составления адекватной математической модели процесса фильтрации подземных вод необходимо иметь более или менее полную информацию о пространственном распределении значений искомых параметров. Определение гидрогеологических параметров опытнофильтрационными методами сопряжено со значительными материальными и временными затратами. В данной работе рассматривается метод идентификации коэффициента фильтрации покровного слоя расчетным путем по некоторым известным значениям уровня грунтовых вод (УГВ).
Задача идентификации коэффициента фильтрации в неоднородной пористой среде сводится к решению коэффициентной обратной задачи для уравнения Бус-синеска. Для обеспечения единственности решения должны быть заданы значения УГВ и искомой функции в некотором дискретном множестве точек, полученные наблюдением и/или экспериментом.
Установившееся движение грунтовых вод в пористой среде, ограниченной сверху поверхностью земли, а снизу - водоупором, описывается уравнением Бус-синеска [1]:
Т —1_ —ГТ —1 = /(х,у), (х,у) е Б (1)
дх I йх) ду ^ ду ) ^
с граничным условием
Т — = а+рЯ, (х,у) еГ = дБ, (2)
дп
где
Т = к(Н _ Ь) _ водопроводимость, м2/сут.; (3)
И = И(х, у) _ УГВ, м; к = к(х, у) _ коэффициент фильтрации покровного слоя подлежащий идентификации, м/сут.; Ь = Ь(х, у) _ поверхность водоупора подстилающего покровный слой, м; /(х, у) _ функция инфильтрации и испарения, м/сут.;
а = а(х, у) и р = Р(х, у) _ заданные функции; Б _ область фильтрации в плане, Г _ ее граница.
В данной задаче, кроме условия (2), задаются, так называемые, внутренние граничные условия
о которых говорилось выше.
Задача заключается в определении функции Т(х, у) из уравнения (1) при соблюдении условий (2) _ (5). Поскольку значения УГВ задаются с определенной погрешностью и в недостаточном объеме, то задача нахождения коэффициента уравнения (1) является некорректной, поэтому для ее решения применяем метод регуляризации А.Н. Тихонова [2]. Задача сводится к нахождению функции Т(х, у), сообщающей в области Б минимум функционалу [3]
где 5Т _ вариация функции Т(х, у); у _параметр регуляризации; Н(Т) _ расчетные значения УГВ, которые находятся как решение задачи (1), (2).
Задача (1), (2) решается методом конечных элементов [4]. Область фильтрации Б разбивается на треугольные элементы таким образом, чтобы точки, в которых заданы экспериментальные значения И? и К]■, совпали с вершинами элементов. В этих точках используются условия (4) и (5), а в остальных узлах начальные приближения функций Н(х, у) и К(х, у) должны удовлетворять условиям
шшН < И < тахН и ттКэ < К < тах Кгэ соответственно. Начальное рас-
I г ] у
пределение функции Т(х, у) вычисляется по формуле (3).
В элементе (е) с вершинами г,], к функция Н(х, у) выражается формулой [5]
H(x, yt) = Hi , i = 1,2,..., p,
K(xj,yj) = Кэ , j = 1,2,...,q,
(4)
(5)
p 2 q
Ф(т) = £[ Hi (T) - Hi ] + £ (Tj - tj )2 +y|st| 2,
(6)
j=1
H(e) (x, y) = HiNie) (x,y) + HjNf (x, y) + HkNe (x,y),
(7)
где
Hs = H (xs, ys ), Ns ( x, y) = (as + bsx + csy)l2A, s = h j, k, ai = xjxk -xkyj, bi = yj -yk, ci = xk -xj,
aj = xkxi - xiyk, bj = yk - yi, cj = xi - xk,
ak = xiyj - xjy, bk = yi - yj, ck = xj - xi,
1 xi yi
2A = 1 xj yj .
1 xk yk
Суммируя равенство (7) по всем элементам (е), получаем формулу
П
(8)
j=1
где п _ число всех узлов сетки.
Применяем к задаче (1), (2) принцип Галеркина:
-Ц Ы,(х У)
и
—( Т —']+—Г Т — ] + / дх ^ дх ) ду ^ дх )
ёхёу + Г Ы, (х, у) ^ Т дП _ а _ рн jd.s• = 0,
] = 1,2,..., и.
Используя в двойном интеграле формулу Грина, приходим к системе уравнений
Я
и
дН дЫ, дН дЫ,
дх дх ду ду
ёхёу _Ц Ыjfdxdy _| Ы,Н вds = 0, і = 1,2,..., п.
и
После подстановки вместо Н(х, у) ее разложения (8) получаем систему линейных алгебраических уравнений относительно УГВ
£ ачн} =ь, ]=і
] = 1,2,..., п,
(9)
где
а] =11Т
и
дЫ, дЫ] дЫ, дЫ] ]
дх дх ду ду
Ь, = И Ы, fdxdy +1Ы, аds .
dxdy _ 1Ы, Ы ] р ds,
(10)
и
Систему (9) решаем одним из точных или приближенных методов. Подставляя полученные значения Н(х, у) в формулу (3), решаем систему (9) при новом значении Т(х, у) и т.д. до выполнения условия
тах
1<і <п
н(у) _ н(у-1)
< є,
(11)
где V - номер итерации; е > 0 - заданное малое число.
Для определения поля функции Т(х, у) мы наряду с количественной информацией (условия (4) и (5)) используем также качественную информацию об искомой функции, т. е. функционал Ф(Т) требует, чтобы функция Т(х, у) была гладкой, что соответствует физической природе водопроводимости.
Теперь займемся минимизацией функционала (6). При каждом наборе значений функции Т(х, у) получаем вполне определенные значения функции Н(х, у), т. е. имеем оператор Н(Т), определенный алгоритмически по формулам (1), (2), (9), (10). Этот оператор в общем случае является нелинейным. Линеаризуем его следующим образом
Н(Т) = Н(Т) + £(Т, _Т,)-|Н + Я2(АТ);
,=1 ■«
(12)
где Тв -значение функции Т в точке (хл у), полученное в предыдущей итерации; Я2 (АТ) - остаточный член разложения. Подставляя (12) в (6) и используя необходимое условие минимума функции многих переменных, получаем систему линейных алгебраических уравнений
дФ(Т) = £
дТк
і=1
дН,
Н,+£ (Т _ Т,) -Т- _ Н,
^ + цТ _ Тк) + у(Тк _ Тк) = 0,
дТк
к = 1,2,..., п,
или
с коэффициентами
VCksTs - dk , k - I,2,...,
s-1
(13)
Cks -V dH^, k * s ,
Й dTs dTk
Ckk - £
и с правыми частями
-\dTk j
,dH,
+ у + ц
dk -V Hf - Hi + Y^-Ts
k Z~1 1 dTs s
i-1V s-1 s j
dH,
Ж
■ + M-Tk3 +У Tk.
Здесь
Ц-'
|1, если kэ задано,
|0, если kэ не задано.
Матрицы систем (9) и (13) являются симметричными и имеют диагональное преобладание, и они легко решаются методом Гаусса. Производная дН / дТ является операторной, т. е.
дН дТ
3T, v J j
i, j - 1, n ,
а частные производные аппроксимируются разностными отношениями
дНг _ АН,- = Н,(у+1) - Н,(у) дТ, ~ АТ, ~ Т(v+1) - НМ ,
где V - номер итерации.
Вычислительная процедура осуществляется в следующем порядке. Используя начальные значения Н и Т в качестве нулевого приближения, решается задача (1), (2) и определяется первое приближение Н(1). Затем, придавая приращение Ак функции к(х, у), находим второе приближение Н(2). Это дает возможность приближенно определить производную дН / дТ и решить систему (13) при некотором значении параметра регуляризации у. Итерация по V проводится до установления фильтрационного процесса. Если при этом полученные значения УГВ в
пределах ошибок не совпадут с данными экспериментальными значениями Нэ, то итерация проводится по параметру у, выбор которого может быть осуществлен методом невязок [2].
Алгоритм идентификации коэффициента фильтрации проверен на решении ряда тестовых задач. Рассмотрим одну из них. В центре круговой области с радиусом I = 3000 м расположен постоянный источник с заданным дебитом, а также происходит площадное питание грунтовых вод с заданной интенсивностью. Фильтрация происходит от центра области к ее границе, где задано краевое условие (2). УГВ меняется от 370 м (в центре области) до 350 м (на границе), а коэф-
Точные и приближенные значения функции к(х,
Точные
Приближенные значения коэффициента фильтрации
значения коэффи- p = 39 p = 15 p = 10
циента фильтра- ции q = 15 q = 10 q = 5 q = 15 q = 10 q = 5 q = 15 q = 10 q = 5
11,63 12,85 10,62 11,63 11,99 10,62 11,63 11,99 10,62 11,63
15,36 16,74 17,05 16,55 16,14 17,05 16,55 17,50 17,95 18,61
22,62 22,62 25,04 25,85 22,62 25,84 23,80 22,62 26,45 23,80
35,04 34,01 36,95 37,04 34,01 36,95 31,01 34,01 36,95 31,01
50,00 50,00 55,98 51,60 55,00 54,90 51,60 55,87 56,97 57,95
61,07 62,65 63,88 60,20 62,65 63,88 60,20 62,65 63,88 60,20
64,96 66,98 64,96 60,55 72,12 64,96 74,85 73,56 64,96 76,98
77,38 75,95 74,02 70,37 75,95 74,02 75,99 75,22 74,02 79,94
84,64 84,64 89,00 80,55 84,64 89,00 80,55 84,64 89,05 80,55
88,36 89,70 85,05 82,31 89,70 85,05 81,02 88,01 89,94 81,02
Относи-
тельные погреш- ности 10,5% 12% 14,3% 11% 14,2% 15,2% 14% 16,9% 21,2%
фициент фильтрации - от 10 м/сут до 90 м/сут. Функции, фигурирующие в задаче (1), (2), задаются по формулам:
H(х,y) - ij3502 + e(l2 - x2 -y2) , k(x,y) - 10[5 + 4sinу (x + y)], b(x,y) - 0,
f (x, y) - 2e
20п n
k (x, y) + — (x + y)cos у (x + y)
+ QH, P(x, y) -
e l k(x, y) + a(x, y) H (x, y)
где 0>(х, у) и а(х, у) - произвольные функции, е = 16-10""4.
Задавая в определенных точках области экспериментальные (точные) значения УГВ и коэффициента фильтрации, по описанной методике идентифицируем функцию к(х, у). В табл. 1 приведены точные и приближенные значения функции к(х, у), полученные при задании экспериментальных значений этой функции в q = 15,10,5 точках. При этом значения УГВ задаются в р = 39,15,10 точках, а общее количество узлов сетки п = 39.
Список литературы
1 Полубаринова-Кочина, П.Я. Математические методы в вопросах орошения / П.Я. Полубаринова-Кочина, В.Г. Пряжинская, В.Н. Эмих. - М.: Наука, 1969. -414 с.
2 Тихонов, А.Н. Методы решения некорректных задач / А.Н. Тихонов, В.Я. Арсенин. - М.: Наука, 1974. - 233 с.
3 Мурзакматов, М.У. Об идентификации коэффициента фильтрации в неоднородном водоносном горизонте / М.У. Мурзакматов, Ж.М. Мамыров // Вестник Иссык-Кульского университета. - №3. - 1999. - С. 73-77.
4 Джаныбеков, Ч. Математическое моделирование движения грунтовых вод в многослойных средах / Ч. Джаныбеков. - Фрунзе: Илим, 1982. - 280 с.
5 Сегерлинд, Л. Применение метода конечных элементов / Л. Сегерлинд. -М.: Мир, 1979. - 392 с.
Identification of the Filtration Coefficient of Soil Ground by Regularization Method
M.U. Murzakmatov, B.A. Baibolotov
Issyk-Kul State University (Kyrgyz Republic)
Key words and phrases: filtration coefficient; identification; regularization; finite elements method; subsoil waters.
Abstract: Application of A.N. Tikhonov’s regularization method to identification of filtration coefficient in two-dimension Bussinesk’s equation is considered.
Identifizierung des Koeffizienten der Filtrierung der Boden durch die Methode der Regulation
Zusammenfassung: Es wird die Anwendung der Methode der Regulation von A.N. Tihonov zur Identifizierung des Koeffizienten der Filtrierung in der zweidimensionalen Bussinesk-Gleichung betrachtet.
Identification du coefficient de la filtration des sols par la methode
de la regularisation
Resume: Est etudiee l’emploi de la methode de la regularisation de A.N. Tihonov pour l’identification du coefficient de la filtration dans une equation bidimensionnelle de Bussinesk.
