УДК 519.632.4:622.323
© Н.Д. Морозкин, Г.С. Бикбулатова, 1998
Н.Д. Морозкин, Г.С. Бикбулатова (Башкирский государственный университет)
N.D.Morozkin, G.S.Bikbulatova (Bash State University)
Расчет пластового давления с помощью метода конечных элементов
Application of finite elements' method for calculation of layer's pressure
For the solution of problem of distribution of layer's pressure finite diffrence method is mostly used. In this article the solution of this problem is found by means of finite elements' method. Comparing to the finite difference method, finite elements' method makes it possible to work with a more large net, at the natural account of 6-function, entering the equation and modelling the activity of forcing and extracting wells. The results of the calculating experiment which confirms the effectiveness of the proposed appoach are fuggested.
Г г
ри решении двумерных задач расчета пластового давления в нефтяных месторождениях обычно используют метод конечных разностей. В данной статье для их решения применяется метод конечных элементов (МКЭ). По сравнению с методом конечных разностей он позволяет работать с более крупной сеткой (при той же точности) за счет естественного учета 6-функций, входящих в уравнение и моделирующих работу нагнетательных и добывающих скважин. Приведены результаты вычислительного эксперимента, подтверждающие эффективность указанного подхода.
Распределение пластового давления р(х, у) в некоторой плоской области й£Л2 описывается уравнением [1]
-с^[а(х,у^гаСр(х,у)] =
т
= 2 Яр^х.^ у-у) (1)
где а у) = к Н (х у) F (х y),
а(х,у)>С0>О, к - тензор абсолютной проницаемости; Н(х,у) - нормализованная к единице переменная толщина; F(x,y) - доля нефти в потоке или функция Баклея-Леверетта; 6(х-х., у-у.) -6-функция Дирака; д. - константы, характеризующие мощность скважин (для добывающих скважин д<0, для нагнетательных д. >0); (х, у.), I = 1,2,..т - координаты скважин; т - число скважин.
Сумма в правой части уравнения (1) характеризует работу нагнетательных и добывающих скважин, которые представлены в виде точечных источников, что на практике вполне приемлемо, поскольку диаметр скважин достаточно мал по сравнению с размерами месторождения.
Предполагается, что на границе 5 области й давление задано
Р(х,У) 18 = g(x,У). (2)
Кроме того, для корректности последующих математических выкладок будем считать, что функции а(х,у) и g(x,y) непрерывны и имеют кусочно-непрерывные производные соответственно в области й и на границе 5.
В научной литературе для поиска решения уравнения (1) при граничном условии (2), как правило, используется метод конечных разностей [1]. При этом точечные источники перемещаются в соответствующие узлы сетки. В результате для получения приемлемой на практике точности решения необходимо использовать мелкую сетку, что приводит к трудоемкому алгоритму в отношении как времени счета, так и памяти компьютера.
Следует также подчеркнуть, что задача разработки эффективного алгоритма решения системы (1), (2) особенно актуальна при расчете нефтенасыщенности
нефтяного месторождения. В этом случае строится разностная схема во времени и на каждом временном шаге решается задача расчета пластового давления
р(х,у).
Нами для нахождения решения системы (1), (2), т.е. для расчета давления р(х,у) в нефтяном месторождении, используется МКЭ. Этот метод при негладких входных данных, как правило, сходится быстрее, чем метод конечных разностей, а иногда обладает оптимальной скоростью сходимости [2].
Поскольку правая часть уравнения (1) разрывная функция, решение задачи (1), (2) будем принимать в обобщенном смысле, т.е. под ее решением понимаем функцию и(х,у) £ ^21(^), удовлетворяющую тождеству
гг ( + др
/ у)(--+--) ДО =
й дх дх ду ду
= II qp(x - X y -
Qi=1
(3)
E W24Q),
для любых ^(х,у)
Мх,у) |8 = 0.
Можно показать, что условия на данные задачи (1), (2) обеспечивают ее однозначную разрешимость в ^21(й) [3].
Согласно методу конечных элементов [4] задачу поиска обобщенного решения уравнения (1) с краевым условием (2)
28 11/1998
Таблица 1
Номер скважины x, м У, м k, мкм2 H, м Q,т/сут Sh
1 143239 131630 449,6 14,4 0,399 -169,4
2 143418 132364 364,1 6,8 0,431 -103,5
3 142613 132071 713,3 12,4 0,388 0,0
4 141906 130981 683,6 20,0 0,425 0,0
5 141893 131647 102,3 4,8 0,470 0,0
6 141947 132490 40,5 14,0 0,529 0,0
7 141023 132185 409,2 16,0 0,781 0,0
8 140912 130425 336,2 8,0 0,781 0,0
9 140878 131105 428,4 12,0 0,781 0,0
10 140901 131647 537,1 10,4 0,789 256,0
11 143305 130130 690,9 9,4 0,791 0,0
12 143613 131052 1000,0 5,6 0,512 -84,2
13 142622 130544 223,2 6,4 0,563 0,0
14 141804 130139 420,4 14,8 0,779 0,0
15 143780 130127 818,4 6,6 0,792 117,5
16 143879 132428 759,0 8,4 0,789 185,9
17 142642 132441 432,2 13,0 0,774 488,0
18 141384 130260 203,1 11,6 0,768 84,4
19 141610 132511 531,0 20,4 0,789 0,0
20 144134 130448 471,6 8,6 0,782 176,2
21 144284 131573 811,7 10,8 0,791 352,4
22 142897 130500 96,1 14,8 0,504 0,0
23 143029 131133 546,2 8,2 0,409 -40,3
24 142965 132020 166,9 8,0 0,519 -145,8
25 142696 130884 549,2 6,0 0,498 -76,9
26 142209 130494 781,6 12,6 0,481 -125,2
27 142169 131334 187,1 9,4 0,530 -6,1
28 142245 132010 206,1 9,2 0,452 -185,6
29 143338 131984 381,4 14,8 0,449 -266,9
30 143024 132370 597,9 12,4 0,529 -59,4
31 142554 132412 126,0 13,6 0,524 -86,3
32 143885 132036 650,8 10,4 0,459 -98,1
33 142758 131526 966,0 12,0 0,513 -213,5
34 143606 131501 346,8 8,0 0,457 -120,9
35 142323 130896 198,4 10,0 0,490 -79,1
36 143567 130575 358,0 5,2 0,541 -67,4
37 143240 130857 212,5 7,4 0,479 -77,2
Примечание. Q - продуктивность скважин, SH - нефтенасыщенность. Константы q; вычислялись по формуле
сведем к поиску решения вариационнои задачи о минимуме функционала
ДР) = IНх,уЖ^)2 +(дР)2] -
й дх ду
- 22 P(x,y)qp(x - x.t,y - y.t)}dQ i=1
(4)
на множестве функций из р(х,у), удовлетворяющих условию (2). Вычисляя первую вариацию, получим условие минимума этого функционала 0Т П ! Мдр + др
о/ = / Мх у)[--+--] -
й дх дх ду ду
Таблица 2
кв 0,000 0,006 0,015 0,027 0,041 0,054 0,065
кн 0,892 0,726 0,581 0,461 0,365 0,285 0,220
SH 0,320 0,343 0,366 0,389 0,412 0,435 0,458
кв 0,073 0,082 0,092 0,105 0,117 0,132 0,147
кн 0,165 0,117 0,080 0,054 0,035 0,022 0,013
SH 0,481 0,504 0,527 0,550 0,573 0,596 0,619
кв 0,164 0,182 0,199 0,217 0,235 0,252 0,270
кн 0,008 0,005 0,004 0,003 0,002 0,001 0,000
SH 0,642 0,665 0,688 0,711 0,734 0,757 0,780
/ Г \
-2Ф(х,уШх - x„y - y,)}dQ = 0, (5)
i=1
где ф - произвольная функция из W2l(ß), обращающаяся на границе области Q в нуль.
Будем искать приближенное решение уравнения (5). С этой целью область Q разобьем на конечное число треугольных подобластей (элементов). Вершины треугольников (xj, У;) назовем узлами сетки. В пределах каждого треугольника A(n) приближенное решение р (n) (x,y) будем искать в виде линейной функции ~(n) (x, y) = ax+by+c, (6)
а во всей области Q - в виде суммы функций р(n) (x, у).
Приняв в качестве неизвестных величин параметры P;, i=1, N в узлах сетки (N - число узлов), константы a, b, c для треугольного элемента A(n), определяемого вершинами i, j, k, получим из следующей системы линейных уравнений f ax.+by.+c=P.
{ ax+by+c=Py (7)
ахк+ЬУк+с=Рк
Решив систему (7) и подставив найденные константы a, b, c в выражение (6), приближенное решение уравнения (5) во всей области Q можем записать в виде
[4, 5]
N
p(x, у) = 2 № (х, у), (8)
i=i
где pj, i = 1, N - давления в узлах сетки; Ф^^З^а^у+у;), aj=yj-yk,
ßJ=xk-xj, YJ=xjУk-xkУj, 3(n) - площадь треугольного элемента A(n).
Подставив вместо p(x,y) в уравнение (5) функцию р(х,у) и приняв параметр ф последовательно равным ф1, Ф2,.....,ф^, получим систему для определения неизвестных узловых значений P|,
l = О [4,5], 2 [c(n)p(n) - f 2Фг(x,y)qßx
nEL An s=1
x (x - xs,у - ys)dxdy] = 0, г = 1, N, (9)
где символ n E L означает, что в уравнении с номером г входят элементы, примыкающие к l-узлу; P|(n) = (Pj, Pj, P^)T;
11/1998 29
Р., Р^, Рк - давление в вершинах треугольного элемента Д(п); Т - знак транспортирования; С|(п) = (с|.(п), с^(п), С|к(п));
1
clq
7(ala q + PiPq ) / a(x, y)dxdy;
(2 5(п))2
Система уравнений (9) представляет собой систему линейных алгебраических уравнений для определения давления в узлах сетки. Окончательный вид этой системы записывается с учетом краевого условия (2).
Приведем результаты вычислительного эксперимента для решения следующей задачи: определить пластовое давление в нефтяном месторождении площадью 7,8 км2 с 18 добывающими, 7 нагнетательными и 12 отключенными скважинами.
На границе области задавалось постоянное давление, равное 26,3 МПа, т.е. p(x,y) 18 = 26,3 МПа.
Исходная информация по каждой скважине приведена в табл. 1.
Константы д. вычислялись по формуле
д.^/рД., ¡=1,2.....т, (10)
где Qi - продуктивность ¡-й скважины; р - плотность; Д. - площадь поперечного сечения ¡-й скважины. Функция Баклея-Леверетта записывалась в виде [1]
F(x,y)=[kн(Sн)/^н]+[kв(Sн)/^в], (11)
где kн(Sн), кв(5н) - относительная фазовая проницаемость соответственно для нефти и воды; = 2,59 мПа с, = = 0,42 мПа • с - вязкость соответственно нефти и воды; 5н = 5н(х,у) - нефте-насыщенность в точке (х,у).
Функции кн(5н) и кв(5н) задавались таблично в 21 точке. Значения параметров приведены в табл. 2.
Нефтенасыщенность в узлах сетки пересчитывалась по известным значениям нефтенасыщенности в скважинах по формуле
5: = ^- , ]= 1, 2.....N
j 2 г
(12)
где r - расстояние от j-го узла до г-и скважины. Суммирование проводилось
по скважинам, попадающим в круг некоторого характерного радиуса с центром в j-м узле. Аналогично пересчитывались значения к и H(x,y).
Проницаемости кн и кв в узлах сетки определялись следующим образом. Из условия S <S <
И1 н
<S (1+1), г = 1,2.....21,
н(1+1).....
j = 1,2,..., N выбирался интервал (S^, SИ(1+1)). По приведенным выше значениям к ,S , к определялись интервалы
[kи(Sиl), №(!+!))]
и [kв(Sиl), kв(Sи(1+1))],
в которых параметры к ( S ) и к ( S ) аппроксимировались линеино.
Разбиени
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.