CC BY
50
УДК 517.958:532 ББК 22.181 Ш 96
М.М. Шумафов, Р. Цей
Разработка алгоритма для численного решения обратной задачи теории фильтрации методом модулирующих функций
(Рецензирована)
Аннотация
В статье разработан алгоритм определения фильтрационных и ёмкостных параметров газоносного пласта. Предложенный алгоритм основан на решении обратной задачи теории фильтрации методом модулирующих функций. Также приводится численная реализация решения поставленной задачи.
Ключевые слова: обратная задача, математическое моделирование, идентификация, модулирующие функции, теория фильтрации, фильтрационные и ёмкостные параметры.
M.M. Shumafov, R. Tsey
Development of an algorithm for a numerical solution of an inverse problem of the filtration theory by modulating function method
Abstract
The paper suggests the algorithm for determining filtrational and capacity parameters of gas reservoir. This algorithm is based on the solution of an inverse problem of the filtration theory by modulating function method. A numerical solution of the problem stated is given as well.
Key words: inverse problem, mathematical simulation, identification, modulating functions, filtration theory, filtrational and capacity parameters.
Введение
В данной работе нами предложен алгоритм (с численной реализацией) определения фильтрационных и ёмкостных параметров газоносного пласта на основе решения обратной задачи теории фильтрации. Для решения этой задачи применяется метод модулирующих функций (далее, М-метод).
Отметим, что идея применения М-метода для решения обратных задач восходит к работам Дж. Лоэба и Г. Кахена (J. Loeb, G. Cahen) [1, 2]. Возможность применения М-метода для решения задач нефтегазовой науки впервые была высказана В.Б. Георгиевским и им были разработаны унифицированные алгоритмы для решения обратных задач подземной гидрогазодинамики [3]. В работах [4,5] приведены некоторые способы программной реализации алгоритмов, предложенных в работе В.Б. Георгиевского [3]. В работе [6] сделана попытка обобщить М-метод на случай любой степени полиномов разложения неизвестных параметров газоносного пласта.
Решение обратной задачи теории фильтрации
В работе [7] рассмотрена следующая обратная задача теории фильтрации.
Дано нелинейное эволюционное уравнение параболического типа, описывающее процесс нестационарной фильтрации реального газа [8]
3 а( х, у) 3 р2 3 + 1 аа хх )у 3 р2 = Ъ(х, у) / \ р
3 х т (р)г(р) 3 х . 3 у )р N )р 3 у . 3г 2( р) \ /
+ с( х, у, г),
(х, у) е О , г е (0,Т), п с я2 с начально-граничными условиями
Р(хУ,0) = Ро(XУ), (хУ)е п ,
Lp| = а (х, у, г), (х, у) е Г , г е [0, Т],
(2)
(3)
где а( X, у) = k (X, у)А( х, у), Ъ( х, у) = 2а т( х, у)А( х, у), с( х, у, г) = 2Q( х, у, г) рат.
Здесь р( х, у, г) - давление в точке пласта с координатами (х,у) в момент времени г, к(х,у) - коэффициент проницаемости пласта, т(х,у) - коэффициент пористости пласта, ^х,у) - эффективная толщина пласта, /и(р) и г(р) - соответственно коэффициенты динамической вязкости и сверхсжимаемости газа при давлении р и пластовой температуре Тпл, а(х,у) - коэффициент газонасыщенности, Q( х, у, г) - объемный расход газа, отнесенный к единице площади пласта в точке (х,у) в момент времени г, приведенный к атмосферному давлению рат и пластовой температуре Тпл. Далее, L - это единичный оператор в случае первой краевой задачи; в случае второй краевой задачи L = 3 / 3 п - производная по внешней нормали к границе Г области О ; в случае третьей краевой задачи Lp = 3р/3п + Q(х,у,г)р (Q(х,у,г), а (х,у,г) - известные функции).
Отметим, что оператор L в условии (3) может иметь и более сложную структуру.
Требуется найти коэффициенты а( х, у), Ь(х, у) в уравнении (1) при условии что функции т (р), 2(р), с(х, у, г) известны.
Введем в рассмотрение функцию от давления
2X4
У (р): ]
0 т (X) ^ (X)
тогда уравнение (1) можно переписать так [7]:
_3_ 3 х
а( х, у)
3х
+
_3_ 3 у
а( х, у) 3'У~ 3у
= Ъ(х,у)^ +у 2, 3г
(4)
(5)
где У ! = р / z(p), У 2 = 2рат^X, У, г) .
Для нахождения коэффициентов а(х, у) и Ъ(х, у) в (5) применяется М-метод. Разлагая коэффициенты а(х, у) и Ъ(х, у) по формуле Тейлора
k- туШ +
a(х, у) = XX а
k = 0 т= 0
Ъ(х,у) =Хр X Ъкпхк~тут + ...
(6)
(7)
k = 0 п= 0
после элементарных преобразований получим следующую систему алгебраических уравнений
X X (а ЙЧ, + Ь “Ъ„) = 7 “ и = °.1.-...«.- 1),
k = 0 т= 0
или в векторно-матричной форме
Ж'0 = и . (8)
Здесь Ж - матрица размера пе х пе з элементами которой являются известные коэффициенты « (кл, Р(кл, а 0 - ", -мерный вектор неизвестных коэффициентов , Ъкт .
Нетрудно установить, что при к= 0,1 ,...,прг1 т- 0,1,...,А: число неизвестных
{'акт) и равно, соответственно,
(пр + IXй, + 2)
п„ =
и , а общее число
с помощью индекса /
неизвестных равно пе - (пр + 1 )(пр + 2).
Вместо двойных индексов кт мы иногда будем использовать одинарный индекс /
(схема нумерации с помощью индекса / изображена на рис. 1). Таким образом, индекс
будет обозначать номер столбца, а / - номер п , г , ,
■’ ' ^ 1 рис | Схема нумерации коэффициентов
строки матрицы Ж.
С помощью индексов /, j элементы матрицы Ж и векторов 0 и У запишутся в виде:
^ ■= {а,(Л,Р,(Л}, 0 = {аМТ, У = {/«},
(7= 1,2,...,Ив; j= 0,1,...,Ив - 1).
Здесь знак Т означает операцию транспонирования.
Числа й г(;), Р -л, / (;) определяются, соответственно, следующим образом:
= «£> = ( /А<)ШуЛх‘'У"/;(х)* (к~ т)х‘-'у"Г,(х)\+
V
+ Ш*)[^ ^7ГОО+
(9)
РУ ■= № -- I¥г/]ш](уу~тут/,](^г
V
г°’= |»',Л«Л(г)Л(0Л'
* 5
(10)
(П)
где К = {[Хі,х2]х (к= 0,1,.~,пр; т=0,1,.~,к;
7= 0,1,..., пе- 1).
Для однозначного определения пе неизвестных аы з Ъкт з необходимо иметь >>,, -уравнений, причем детерминант матрицы Ж линейной системы (8) должен быть отличен от нуля.
Для этого модулирующие функции /, 00,./, О’Х ./, (0 выберем следующим образом: /Дх) = х;(х- х1)2(х2 - х)2, хє [х15х2],
ї](у)= у](у- Уі)2(У2- у)\ ус [Уі,У21, /.(0= *4*- 0, ^ ІАЛ],
(12)
(7= 0,1,...,ив - 1).
Из формул (9)- (11) видно, что вопрос машинной реализации поставленной задачи сводится к вопросу численного интегрирования тройных интегралов.
Алгоритм для численного решения обратной задачи теории фильтрации
Из сказанного выше получаем следующий алгоритм численной реализации для решения обратной задачи теории фильтрации.
Шаг 1. Определение числа неизвестных, а тем самым и порядка матрицы Ж
п, = (пр + 1)(пр + 2),
где пр - максимальный порядок степеней разложения искомых функций а(х, у) и Ъ(х, у) в разложениях (6), (7).
Шаг 2. Генерирование модулирующих функций fj(х), fj (у), fj(г),
удовлетворяющих условиям (12). Для этого удобно представить модулирующие функции следующим образом:
х 4 у 2 г
а х (у) = у1 X а у , (г) = IіX а г , (1 = 0,1>->п -1),
V
V = 0
V = 0
V = 0
х у Ґ
где соответственно равны
х х
а 0 = х1 х2, а 1 = - 2(х1 + х2)х1 х2
а 2 = (х1 + х2)2 + 2х1 х2, а 3 = - 2(х1 + х2), а 4 = 1,
у у
а 0 = у"у2 , а 1 =- 2(у + у2)у1 у2’
у у у
а 2 = Сух + у 2) + 2у1 у2, а 3 =-2(у1 + у2), а 4 = 1,
г г г
а 0 = - іхі2, а х = (гх + г2)3 а 2 = -1-
Шаг 3. Получение и обработка входных данных. Для решения задачи сначала необходимо иметь экспериментальные данные о значениях р( х, у, г), Q( х, у, г) в узлах кубической решетки, а также необходимо вычислить /л(р) и z(p).
При отсутствии данных во всех узлах решетки можно интерполировать функции р(х, у, г), Q(х, у, г). Для этого можно применить известные методы интерполяции (например, интерполяционные многочлены Лагранжа, Ньютона и др.) [9-12].
В соотношениях (9)-(11) присутствуют У = У (р), у х = у х(р) , у 2 = у 2(х, у, г) .
Сначала вычислим У (р) из (4), взяв в качестве z(p) квадратичную функцию:
z(P) = а0Р2 - Ь0Р + С0, (13)
где а0 = 0.18348• 10-5, Ь0 = 0.50654• 10-3, с0 = 0.9911. В качестве т (Р) возьмем ее усредненное значение, т.е. положим т (р) = т ср . Так как 4а0с0 - Ь02 > 0, то
2 pdp
У (Р):
а0 Р - Ь0 Р + С0
0 ^0
1
а
^п(а0Р2 - Ь0Р + С0) +
2^0
ат^
2а0Р - Ь0
■10Н ^0
,2 7.2
+ с,
где с = - —
т 0
— 1п с0 + а
2^0
ат^-
I------2----і-----------------------_
д/4а0С0 - Ь0 л/4а0С0 - ь
Окончательно имеем:
У (Р) =
1п(а0Р2 - Ь0Р + С0)
2^0
2а0 Р - Ь0 атсг^~і= - ат^-
I-----2-----2 6 I------2-----2
'У4а0С0 - Ь0 л/4а0С0 - Ь(
Р
0
і
+
с
0
Выбор функции р) мотивируется тем обстоятельством, что для глубинных газовых залежей коэффициент сверхсжимаемости хорошо аппроксимируется квадратичной функцией вида (13). Отметим, что квадратичная парабола дает хорошее приближение одной из эмпирических кривых Д. Брауна [13-15] (см. рис. 2).
Рис. 2. Зависимость коэффициента сверхсжимаемости от приведенных давления и температуры для природных газов [13]
Шаг 4. Заполнение значений элементов матрицы Ж и вектора У (см. (8)), а именно значений а 1^, Ь ы, 1(1) (см. (9)-(11)), представляющих собой тройные интегралы.
Вычисление тройных интегралов в соотношениях (9)-(11) сводится к вычислению повторных интегралов с применением известных формул численного интегрирования (формулы Ньютона-Котеса и др.) [9-12].
Шаг 5. Нахождение сопгі(Ж) - числа обусловленности матрицы W:
сопй (Ж)
где 1 тах ,
матрицы
1 тт - максимальное и минимальное, соответственно, собственные значения W. Для преодоления трудностей, связанных с некорректностью плохо
обусловленных систем, необходимо обеспечить «близость» числа обусловленности со^(Ж) матрицы Ж к 1 [16]. Это достигается варьированием модулирующих функций или усреднением исходных данных.
Шаг 6. Решение системы линейных алгебраических уравнений (8).
Подходящим образом выбирая модулирующие функции, можно обеспечить выполнение условия невырожденности матрицы Ж: det Ж Ф 0.
Для решения системы (8) можно применить, например, хорошо известный метод исключения Гаусса с выбором главного элемента ([9-12]).
В результате решения системы (8) найдем вектор неизвестных 0 = {а^, Ът}, первые п, /2 элементов которого - коэффициенты ат, а последующие п, /2 элементов -коэффициенты Ът .
Шаг 7. Подстановка найденных значений акт, Ът в разложения (6), (7) искомых функций а(х, у) и Ъ(х, у) соответственно.
Заключение
В статье приведен алгоритм численного решения обратной задачи теории фильтрации, а именно задачи определения фильтрационно-ёмкостных параметров газоносного пласта. Эти параметры являются коэффициентами дифференциального уравнения, описывающего процесс неустановившейся фильтрации газа в неоднородной по коллекторским свойствам пористой среде.
Алгоритм основан на решении обратной задачи теории фильтрации М-методом.
Отметим, что одним из главных достоинств М-метода (наряду с эффективностью) является его простота для практического применения.
Примечания:
References:
1. Loeb J., Cahen G. Extraction, a partik des enregistrements de mesures, des parametres dynamiques d um system // Automatisme. 1963. №12. P.17-28.
2. Loeb J., Cahen G. More about process identification // Trans. on Automatic Control. 1965. P. 359-361.
3. Георгиевский В.Б. Унифицированные алгоритмы для определения фильтрационных параметров: cправочник. Киев, 1971. 328 с.
4. Лукнер Л., Шестаков В.М. Моделирование геофильтрации. М.,1976. 407 с.
5. Трофимов В.В., Батищева Г.А. Реализация на ЭВМ унифицированных алгоритмов В.Б. Г еоргиевского // Сборник научных трудов ЮжНИИгидротехники и мелиорации. 1976. Вып. 9. С. 111-114.
6. Юдин А.И., Юдина О.К. Расчет фильтрационноёмкостных параметров по промысловым данным эксплуатации газового месторождения // Термодинамика кооперативных процессов в гетерогенных средах. Тюмень, 1985. С. 80-85.
7. Шумафов М.М., Цей Р. Идентификация параметров газоносного пласта на основе решения обратной задачи теории фильтрации //
1. Loeb J., Cahen G. Extraction, a partik des enregistrements de mesures, des parametres dynamiques d um system // Automatisme. 1963. No.12. P.17-28.
2. Loeb J., Cahen G. More about process identification // Trans. on Automatic Control. 1965. P. 359-361.
3. Georgievsky V.B. Unified algorithms for definition of filtrational parameters. Directory. -Kiev, 1971. 328 pp.
4. Lukner L., Shestakov V.M. Modelling of a geofiltration. M., 1976. 407 pp.
5. Trofimov V.V., Batishcheva G.A. Realization with the help of computer of the V.B. Georgievsky unified algorithms. Proc. Southern Research Hydrotechnical and Melioration Institute, 1976, Issue 9. P. 111-114.
6. Yudin A.I., Yudina O.K. Calculation of filtrational-capacity parameters from trade data on exploitation of a gas deposit // Thermodynamics of cooperative processes in heterogeneous environments. Tyumen, 1985. P. 80-85.
7. Shumafov M.M., Tsey R. Determination of gas reservoir parameters on the basis of a solution of an inverse problem of the filtration theory // The
Вестник Адыгейского государственного университета. Сер. Естественно-
математические и технические науки. Майкоп, 2009. Вып. 1 (43). С.33-42.
8. Закиров С.Н., Лапук Б.Б. Проектирование и разработка газовых месторождений. М., 1974. С.39.
9. Бахвалов Н.С., Жидков Н.П., Кобельков Г.М. Численные методы. М., 1987. 600 с.
10. Калиткин Н.Н. Численные методы. М., 1978. 512 с.
11. Копченова Н.В., Марон И.А. Вычислительная математика в примерах и задачах. М., 1972. 368 с.
12. Самарский А.А., Гулин А.В. Численные методы: учеб. пособие для вузов. М., 1989. 432 с.
13. Басниев К.С., Дмитриев Н.М., Розенберг Г.Д. Нефтегазовая гидромеханика: учеб. пособие для вузов. М.; Ижевск, 2005. С. 406.
14. Басниев К.С., Кочина И.Н., Максимов В.М. Подземная гидромеханика: учеб. для вузов. М., 1993. 416 с.
15. Механика насыщенных пористых сред / В.Н. Николаевский, К.С. Басниев, А.Т. Горбунов, Г.А. Зотов. М., 1970. 339 с.
16. Форсайт Дж., Молер К. Численное решение систем линейных алгебраических уравнений: пер. с англ. М., 1969. 168 с.
Bulletin of the Adyghe State University. Ser. Natural-Mathematical and Technical Sciences. Maikop, 2009. Issue 1(43). P. 33-42.
8. Zakirov S.N., Lapuk V.V. Designing and development of gas deposits. M., 1974. P. 39.
9. Bakhvalov N.S., Zhidkov N.P., Kobelkov G.M. Numerical methods. M., 1987. 600 pp.
10. Kalitkin N.N. Numerical methods. M., 1978. 512 pp.
11. Kopchenova N.V., Maron I.A. Calculus mathematics in examples and in problems. M., 1972. 368 pp.
12. Samarskiy A.A., Gulin A.V. Numerical methods: Manual for higher schools. M., 1989. 432 pp.
13. Basniev K.S., Dmitriev N.M., Rozenberg G.D. Oil-and-gas hydromechanics: Manual for higher schools. M.- Izhevsk. 2005. P. 406.
14. Basniev K.S., Kochina I.N., Maximov V.M. Uunderground's hydromechanics: Textbook for higher schools. M., 1993. 416 pp.
15. Mechanics of the saturated porous environments / V.N. Nikolaevskiy, K.S. Basniev, A.T. Gorbunov, G.A. Zotov. M., 1970. 339 pp.
16. Forsite J., Miller K. Numerical solution of linear systems of algebraic equations (Trans. from Eng.). M, 1969. 168 pp.
