Толпаев В.А. Математические модели двумерной фильтрации в анизотропных, неоднородных и многослойных средах
Подождите немного. Документ загружается.
151
скорости фильтрации
µ
−=ϕ
kP
лишь постоянными множителями отличаются от
приведенного давления P, то поверхность ствола скважины в рассматриваемых
расчетах в то же время принимается и за поверхность с постоянным напором H и
за эквипотенциальную. Между тем для скважин, эксплуатирующих пласты боль-
шой мощности и обладающих большим дебитом, предположение о том, что по-
верхность её ствола эквипотенциальна, может приводить к заметным ошибкам в
фильтрационных расчетах.
Впервые задача об учете градиента потенциала вдоль ствола скважины при
расчете ее дебита для фильтров перфорационных конструкций ставилась в [253] и
наиболее широко для совершенной скважины рассматривалась в [26]. Однако в
[26, 253] вывод расчетного уравнения (схема вывода которого приведена также в
[109]) для переменного вдоль ствола скважины расхода был основан на том, что
напор на внешней стенке фильтра, т.е. на стволе скважины, считался постоянным.
Подчеркнем, что такое допущение о постоянстве напора на поверхности ствола
скважины не совместимо с неравномерным распределением скорости фильтрации
по её высоте.
В этом параграфе автор даёт точное решение уравнений напорной фильтра-
ции в осесимметричной постановке для фильтров, создаваемых способом гравий-
ной набивки [134] и принципиальная схема которых приведена на рис.33. Полу-
ченное решение, в котором учитывается неравномерность распределения приве-
денного давления и скорости притока флюида по высоте фильтра, применено для
точного расчета дебита нефтедобывающей скважины.
4.4.2. Уравнения и граничные условия
Фильтрация флюида в областях 1 и 2 на рис.33 считается линейной, следую-
щей закону Дарси
11
gradv ϕ= ;
µ
−=ϕ
Pk
1
1
(1)
22
gradv ϕ= ;
µ
−=ϕ
Pk
2
2
. (2)
152
Поскольку для несжимаемой жидкости уравнение неразрывности для областей 1 и
2 имеет вид
(
)
0vdiv = , (3)
то после подстановки (1) и (2) в (3) для ϕ
1
и ϕ
2
в цилиндрических координатах по-
лучаем одно и то же уравнение - уравнение Лапласа
2,1i,0
z
r
r
rr
1
2
i
2
i
==
∂
ϕ∂
+
∂
ϕ∂
∂
∂
. (4)
В классической постановке задачи о дебите скважины круговая цилиндриче-
ская поверхность АА
1
В
1
В (рис.33) принимается за эквипотенциальную:
const
C
rr
1
C
=ϕ=ϕ
=
. В такой постановке для давления внутри проницаемого стакана
АА
1
В
1
В (рис.33) по всей мощности пласта считается возможным принять гидро-
статический закон. В результате для дебита скважины получается следующая хо-
рошо известная формула Дюпюи
(
)
µ
−
π
=
Π
c
C1
0
r
R
ln
PPbk2
Q
. (5)
На самом же деле цилиндрическую поверхность АА
1
В
1
В фильтра считать эквипо-
тенциальной нельзя. Действительно, движение флюида вдоль вертикальной оси z
внутри фильтра (область 2 с проницаемостью k
2
) обеспечивается за счет градиента
потенциала и, следовательно, ϕ
1
и ϕ
2
на цилиндрической поверхности AA
1
B
1
B по-
стоянными не являются.
Сформулируем теперь граничные условия для уравнений (4) в точной поста-
новке задачи о дебите скважины. На круговой цилиндрической поверхности пита-
ния с радиусом R скорость фильтрации параллельна подошве и кровле пласта и не
зависит от координаты z. Поэтому
const
Rr
1
=ϕ=ϕ
Π
=
. (6)
Поскольку кровля и подошва пласта непроницаемы, то при
Rrr
с
≤≤
0
zz
bz
1
0z
1
=
∂
ϕ∂
=
∂
ϕ∂
==
. (7)
153
Во второй зоне области фильтрации ось Oz ствола скважины является линией тока
течения. Поэтому
0
r
0r
2
=
∂
ϕ∂
=
. Внутри фильтра подошва пласта тоже непроницаема,
поэтому при
с
rr0 ≤≤
0
z
0z
2
=
∂
ϕ∂
=
. (8)
Скорость фильтрации в сечении А
1
В
1
ствола скважины будем рассматривать как
строго вертикальную и имеющую постоянную величину v
0
. Поэтому при
с
rr0
≤
≤
constv
z
0
bz
2
==
∂
ϕ∂
=
. (9)
Чему равна величина v
0
- этот вопрос пока остается открытым, но если v
0
будет
найдена, то объемный дебит скважины (отнесенный к единице времени) найдется
по формуле
0
2
c
vrQ π= . (10)
Решения двух уравнений Лапласа на границе круговой цилиндрической поверхно-
сти AA
1
B
1
B областей 1 и 2 должны удовлетворять условиям непрерывности давле-
ния и нормальной к AA
1
B
1
B составляющей скорости фильтрации. Эти условия, в
силу формул (1) и (2), будут иметь вид:
сс
rr
2
2
rr
1
1
kk
==
ϕ
=
ϕ
(11)
и
сс
rr
2
rr
1
rr
==
∂
ϕ∂
=
∂
ϕ∂
. (12)
После того как сформулированные краевые задачи для уравнений Лапласа (4) бу-
дут решены, можно будет подсчитать среднее значение 〈ϕ
2
〉 потенциала ϕ
2
в круго-
вом сечении A
1
B
1
ствола скважины:
()
()
()
∫∫∫
⋅⋅ϕ=Θ⋅⋅ϕ
π
=ϕ
c
11
r
0
2
2
c
BA
2
2
с
2
drrb,r
r
2
drdrb,r
r
1
. (13)
15
4
Величина 〈ϕ
2
〉, найденная из (13), будет зависеть от v
0
. С другой стороны, в круго-
вом сечении A
1
B
1
бывает известным приведенное давление P
C
, по которому 〈ϕ
2
〉
будет определяться, в соответствии с (2), по формуле
µ
−=ϕ
C2
2
Pk
. (14)
Поэтому на основании уравнений (13) и (14) через известное значение P
C
удастся
вычислить скорость фильтрации v
0
, а затем по (10) и объемный дебит скважины.
Перейдем теперь к решению сформулированных краевых задач.
4.4.3 Расчет потенциала ϕ
1
(r,z)
Решение уравнения (4) для функции ϕ
1
, удовлетворяющее граничным услови-
ям (6) и (7), находится методом Фурье и имеет вид:
() ()
∑
∞
=
π
=ϕ
0n
n1
b
nz
cosrR
~
z,r
, (15)
где
()
rR
~
n
решение дифференциального уравнения
()
(
)
()
0rR
~
rrR
~
r
dr
d
n
2
nn
=λ−
′
, (16)
а
,...2,1,0n,
b
n
n
=
π
=λ (17)
Если в формулах (15) и (16) выделить отдельно слагаемое для n=0, то, с учетом
граничного условия (6), получим:
() ()
∑
∞
=
Π
π
+
+ϕ=ϕ
1n
n01
b
nz
cosrR
~
r
R
lnCz,r
. (18)
В (18) функции
()
zR
~
n
должны удовлетворять уравнению (16) и, для того чтобы
удовлетворялось (6), граничному условию
0R
~
Rr
n
=
=
. (19)
С помощью подстановки
,
b
nr
r
n
π
=λ=ξ
(
)
K,3,2,1n
=
(20)
155
уравнение (16) приводится к виду
() () ()
0R
~
R
~
1
R
~
nnn
=ξ−ξ
′
ξ
+ξ
′′
, которое представляет
собой уравнение для функций Бесселя [131, 135] мнимого аргумента I
0
(ξ) и K
0
(ξ).
Поэтому общее решение уравнения (16) будет иметь вид [135]
()
π
⋅+
π
⋅=
b
nr
KC
b
nr
ICrR
~
0201n
,
(
)
K,3,2,1n
=
. (21)
Так как при r = R должно выполняться условие (19), то C
1
и C
2
в (21) должны быть
подчинены требованию
0
b
nR
KC
b
nR
IC
0201
=
π
⋅+
π
⋅
, откуда получаем
π
⋅=
b
nR
KAC
01
,
π
⋅−=
b
nR
IAC
02
, где A - произвольная постоянная. Таким образом,
функция ϕ
1
(r,z) в первой зоне области фильтрации будет иметь вид:
() ()
∑
∞
=
Π
π
⋅⋅+
+ϕ=ϕ
1n
nn01
b
nz
cosrwC
r
R
lnCz,r
, (22)
где
()
π
π
π
π
=
b
nR
K;
b
nR
I
b
nr
K;
b
nr
I
rw
00
00
n
, (23)
а I
0
и K
0
— модифицированные функции Бесселя нулевого порядка.
4.4.4 Расчет потенциала ϕ
2
(r,z)
Решение уравнения Лапласа для ϕ
2
(r,z) во второй зоне области фильтрации с
проницаемостью k
2
будем искать в виде
() ()
z,rUr
2
1
z
b2
v
z,r
22
0
2
+
−=ϕ . (24)
После подстановки (24) в уравнение Лапласа (4) для U(r,z) снова получим уравне-
ние Лапласа:
0
z
U
r
U
r
rr
1
2
2
=
∂
∂
+
∂
∂
∂
∂
. (25)
Поскольку, как следует из (24),
z
U
b
zv
z
0
2
∂
∂
+=
∂
ϕ
∂
, то для выполнения условий (8) и
(9) придется потребовать, чтобы
0
z
U
z
U
bz0z
=
∂
∂
=
∂
∂
==
. (26)
156
Решение уравнения (25), удовлетворяющее условиям (26), ищем в виде
() ()
∑
∞
=
+
π
ω+
=
1n
n
c
0
const
b
nz
cosr
r
r
lnDz,rU , (27)
где D
0
и const — произвольные постоянные. После подстановки (27) в (25) для
функций ω
n
(r) получается уравнение 0
dr
d
r
dr
d
r
1
n
2
n
n
=ωλ−
ω
, где
b
n
n
π
=λ . Решение
последнего уравнения, совпадающего с (16), уже известно и имеет вид
()
π
⋅+
π
⋅=ω
b
nr
KC
b
nr
ICr
0201n
. (28)
Произвольные постоянные D
0
, C
1
и C
2
подберем так, чтобы 0
r
0r
2
=
∂
ϕ∂
=
, т.к. ось Oz
является линией тока. Из (24), (27) и (28) следует, что
+−−=
∂
ϕ∂
r
D
b2
rv
r
00
2
π
⋅
π
⋅−
π
⋅⋅
π
+
∑
∞
=
b
nz
cos
b
nr
KC
b
nr
IC
b
n
1n
1211
. (Были использованы известные [131,
135] правила дифференцирования модифицированных функций Бесселя I
0
и K
0
).
Для обращения в нуль
r
2
∂
ϕ∂
при r=0 потребуем, чтобы D
0
=0 и
()
(
)
00KC0IC
1211
=
⋅
−⋅ .
Так как I
1
(0)=0, а K
1
(0)=∞, то C
2
=0, а C
1
— произвольная постоянная. Таким обра-
зом, окончательно для ϕ
2
(r,z) получаем следующее представление:
()
+
−=ϕ
22
0
2
r
2
1
z
b2
v
z,r const
b
nz
cos
b
nr
ID
1n
0n
+
π
⋅
π
⋅
∑
∞
=
,
(29)
где D
n
— произвольные постоянные.
4.4.5 Алгебраизация граничных условий сопряжения
Произвольные постоянные C
0
, C
n
, v
0
, D
n
и const в (22) и (29) подберем так,
чтобы оказались выполненными граничные условия сопряжения (11) и (12). Гра-
ничное условие (11) приводит к равенству
()
=
π
⋅⋅+
+ϕ
∑
∞
=
Π
1n
cnn
c
0
1
b
nz
cosrwC
r
R
lnC
k
1
+
π
⋅
π
⋅+
−=
∑
∞
=
const
b
nz
cos
b
nr
IDr
2
1
z
b2
v
k
1
1n
c
0n
2
c
2
0
2
,
которое удобно переписать в виде:
157
bk2
zv
b
nz
cosAA
2
2
0
1n
n0
=
π
+
∑
∞
=
, (30)
где
22
2
c0
c
0
1
0
k
const
bk4
rv
r
R
lnC
k
1
A −+
+ϕ=
Π
, (31)
()
π
−=
b
nr
I
k
D
rw
k
C
A
c
0
2
n
cn
1
n
n
. (32)
Если функцию
bk2
zv
2
2
0
разложить в ряд Фурье по
π
b
nz
cos
, то величины A
0
и A
n
ста-
нут известными, и поэтому (31) и (32) будут уравнениями относительно C
n
, D
n
, C
0
и const. Недостающие уравнения даёт граничное условие (12). Вычисляя
r
1
∂
ϕ
∂
и
r
2
∂
ϕ∂
, получим:
()
π
⋅⋅
π
⋅+−=
∂
ϕ∂
∑
∞
=
b
nz
cosrW
b
n
C
r
C
r
1n
nn
0
1
, (33)
где
()
π
π
π
−
π
=
b
nR
K;
b
nR
I
b
nr
K;
b
nr
I
rW
00
11
n
(34)
и
∑
∞
=
π
⋅
π
⋅
π
⋅+−=
∂
ϕ∂
1n
1n
0
2
b
nz
cos
b
nr
I
b
n
D
b2
rv
r
. (35)
С учетом формул (33) и (35) граничное условие (12) можно будет переписать в ви-
де:
()
=
π
⋅⋅
π
⋅+−
∑
∞
=
1n
cnn
c
0
b
nz
cosrW
b
n
C
r
C
∑
∞
=
π
⋅
π
⋅
π
⋅+−
1n
c
1n
c0
b
nz
cos
b
nr
I
b
n
D
b2
rv
. Из
этого равенства для коэффициентов C
0
, C
n
и D
n
получаем следующие недостающие
уравнения:
b2
rv
C
2
c0
0
= ;
()
n
cn
c
1
n
D
rW
b
nr
I
C
π
= . (36)
158
Из тригонометрического ряда Фурье (30) для коэффициентов A
0
и A
n
находим зна-
чения:
2
0
0
k6
bv
A =
;
()
n
22
2
0
n
1
nk
bv2
A −
π
= . (37)
Теперь из формул (31), (32), (36) и (37) видим, что
(
)
bk12
b2r3v
r
R
lnC
k
1
k
const
2
22
c0
c
0
12
−
+
+ϕ=
Π
; (38)
(
)
(
)
()
()
π
⋅−
π
⋅π
−⋅⋅
=
1
c
0cn2
c
1cn
22
n
cn10
n
k
b
nr
IrWk
b
nr
Irwn
1rWbkv2
D
; (39)
()
()
()
⋅
π
⋅−
π
⋅π
−⋅
π
⋅
=
1
c
0cn2
c
1cn
22
n
c
110
n
k
b
nr
IrWk
b
nr
Irwn
1
b
nr
Ibkv2
C
. (40)
После подстановок (38), (39), и (40) в формулы (22) и (29) для потенциалов ϕ
1
(r,z) и
ϕ
2
(r,z) получаем следующие выражения:
()
+
+ϕ=ϕ
Π
r
R
ln
b2
rv
z,r
2
c0
1
() ()
() ()
∑
∞
=
π
⋅−λ⋅
π
⋅
π
⋅⋅
π
⋅−
π
1n
c
0cn
c
1cn
2
n
c
1
n
2
0
b
nr
IrW
b
nr
Irwn
b
nz
cosrw
b
nr
I1
bv2
,
(41)
где
1
2
k
k
=λ и
()
+
+ϕλ=ϕ
Π
c
2
c0
2
r
R
ln
b2
rv
z,r
(
)
×
π
+
−+
−
2
0
22
0
22
c0
bv2
r
2
1
z
b2
v
b12
b2r3v
() ()
() ()
∑
∞
=
π
⋅−λ⋅
π
⋅⋅
π
⋅
π
⋅⋅−
×
1n
c
0cn
c
1cn
2
0cn
n
b
nr
IrW
b
nr
Irwn
b
nz
cos
b
nr
IrW1
. (42)
4.4.6 Вычисление дебита скважины
Для завершения расчетов остается вычислить среднее значение 〈ϕ
2
〉 по фор-
муле (13). Для этого предварительно запишем ϕ
2
(r,b). Согласно (42) имеем
()
+
+ϕλ=ϕ
Π
c
2
c0
2
r
R
ln
b2
rv
b,r
(
)
×
π
+−
+
2
0
2
0
22
c0
bv2
b4
rv
b12
b4r3v
159
()
() ()
∑
∞
=
π
⋅−λ⋅
π
⋅⋅
π
⋅
×
1n
c
0cn
c
1cn
2
0cn
b
nr
IrW
b
nr
Irwn
b
nr
IrW
.
(43)
Поскольку
()
=⋅⋅⋅
π
=⋅
π
⋅
∫∫
π
b
nr
0
0
22
2
r
0
0
c
c
dxxIx
n
b
dr
b
nr
Ir
()
π
⋅
π
=⋅
π
π
b
nr
I
n
br
xIx
n
b
c
1
c
b
nr
0
1
22
2
c
, то
()
+
+ϕ
λ
=⋅⋅ϕ
Π
∫
c
2
c0
2
c
r
0
2
r
R
ln
b2
rv
2
r
drrb,r
c
(
)
×
π
+−
+
3
c
2
0
4
c0
22
c
2
c0
rbv2
b16
rv
b24
b4r3rv
()
() ()
∑
∞
=
π
⋅−λ⋅
π
⋅
π
⋅
×
1n
c
0cn
c
1cn
3
c
1cn
b
nr
IrW
b
nr
Irwn
b
nr
IrW
.
(44)
Теперь на основании формулы (13) получаем, что
+
+ϕλ=ϕ
Π
c
2
c0
2
r
R
ln
b2
rv
λ
π
++
R
r
,
b
r
,S
r
bv4
3
bv
b8
rv
cc
c
3
2
00
2
c0
, (45)
где
∑
∞
=
=
λ
1n
3
n
cc
n
A
R
r
,
b
r
,S
; (46)
()
() ()
π
⋅−λ⋅
π
⋅
π
⋅
=
b
nr
IrW
b
nr
Irw
b
nr
IrW
A
c
0cn
c
1cn
c
1cn
n
. (47)
С другой стороны, 〈ϕ
2
〉 вычисляется по формуле (14) через заданное значение дав-
ления P
C
в сечении A
1
B
1
добывающей скважины. Это позволяет из (45) найти
уравнение для оставшейся неизвестной величины v
0
:
()
λ
π
+++
λ=
µ
−
R
r
,
b
r
,S
r
b8
r3
b2
4
1
r
R
ln
b2
rvPPk
cc
3
c
3
3
2
c
2
c
2
c0СП2
. (48)
Вычислив из (48) величину v
0
и подставив ее затем в (10), для объемного дебита
скважины, эксплуатирующей пласт с мощностью b, получим:
(
)
λ⋅
λ
π
+++
⋅µ
⋅
−
π
=
−1
cc
3
c
3
3
2
c
2
c
СП1
R
r
,
b
r
,S
r
b8
r3
b2
4
1
r
R
ln
bPPk2
Q
. (49)
160
Отметим, что в частном случае, когда k
2
=k
1
(т.е. когда λ=1), формула (49) даст
формулу дебита скважины, плоским дном вскрывающей продуктивный пласт с
конечной мощностью.
Из (49) вытекает, что влияние неравномерного распределения потенциала
вдоль фильтра АА
1
В
1
В (рис. 1) можно учесть, если в классическую формулу Дю-
пюи ввести поправочный коэффициент β
β
=
0
QQ , (50)
где
1
cc
3
c
3
3
2
c
2
c
c
R
r
,
b
r
,S
r
b8
r3
b2
4
1
r
R
ln
r
R
ln
−
λ⋅
λ
π
+++
=β
, (51)
а Q
0
— находится из формулы Дюпюи (5).
С помощью формул (46), (47) и (51) были вычислены значения поправочного
коэффициента β для различных соотношений R/r
c
, b/r
c
и k
1
/k
2
. Часть результатов
расчетов представлена в таблице 4.2. Расчеты показали, что если отношение про-
ницаемостей k
1
/k
2
<10
-6
, то практически для всех значений R/r
c
и b/r
c
боковую стен-
ку фильтра АА
1
В
1
В скважины при расчетах дебитов можно считать эквипотенци-
альной поверхностью. При этом, если мощность пласта b/r
c
≤500, то относительная
погрешность в расчетах дебитов не превысит 2,4%. Для отношения проницаемо-
стей k
1
/k
2
>10
-5
заметную роль начинает играть мощность пласта. Для пластов ма-
лой мощности, когда b/r
c
<50 и k
1
/k
2
<10
-4
боковую поверхность фильтра скважины
снова можно считать эквипотенциальной и погрешность расчетов при этом тоже
не превысит 2,4%. Для пластов с мощностью b/r
c
>100 в расчетах дебитов уже нуж-
но будет учитывать поправочный коэффициент β, приведенный в таблице. Для
пластов большой мощности неучёт β может привести в расчетах дебитов к замет-
ным ошибкам. Так, например, для k
1
/k
2
=10
-5
, b/r
c
=500 и R/r
c
=1000 расчет дебита по
формуле Дюпюи приведет к относительной ошибке в 22,7%. Приведенная таблица
дает представление о том, в каких ситуациях расчеты дебитов скважин можно вы-