CC BY
24
УДК 550.837
МАТРИЧНЫЙ МЕТОД РЕШЕНИЯ ЗАДАЧИ О СЕЙСМОЭЛЕКТРИЧЕСКОМ ЭФФЕКТЕ ВТОРОГО РОДА В ГЕОЛОГИЧЕСКОЙ СРЕДЕ БИО
В. П. Губатенко1, И. Г. Московский2
1 Саратовский национальный исследовательский государственный университет имени Н. Г. Чернышевского
E-mail: [email protected]
2 Саратовский государственный технический университет имени Гагарина Ю. А.
E-mail: [email protected]
Предложен матричный метод решения задачи о сейсмоэлектри-ческом эффекте второго рода. Рассмотрено возбуждение этого эффекта в горизонтально-слоистой среде точечным источником механических колебаний. Получено аналитическое решение задачи о сейсмоэлектрическом эффекте второго рода для однородного пространства.
Ключевые слова: сейсмоэлектрический эффект второго рода, среда Био, матричный метод, передаточная матрица.
The Matrix Method for Solving the Problem of the Second Kind Seismoelectrical Effect in a Geological Medium Biot
V. P. Gubatenko, I. G. Moskovskiy
The matrix method for solving the problem of Seismoelectrical effect of the second kind is suggested. The excitement of this effect in a horizontally-layered medium by a point source of mechanical vibrations is considered. An analytical solution of the problem of the second kind Seismoelectrical effect for a homogeneous space is obtained. Key words: seismoelectrical effect of the second kind, medium Biot, matrix method, transfer matrix.
DOI: 10.18500/1819-7663-2016-16-4-241 -247 Введение
Сейсмоэлектрическим эффектом второго рода называется явление возбуждения электромагнитного поля в пористых влагонасыщенных горных породах при распространении в них механических колебаний. Электромагнитное поле сейсмоэлектрического эффекта определяется широким спектром петрофизических характеристик горной породы: пористостью, проницаемостью, свойствами поровой жидкости и др. Данная связь определяет стремление использовать сейсмоэлектрический эффект при решении разнообразных прикладных геологических задач, в частности при разработке новых геофизических методов исследования скважин (ГИС) и зондировании околоскважин-ного пространства. Важной составляющей таких методов является возможность математического моделирования измеряемых при зондировании физических полей. Предлагаемая в настоящей
статье математическая модель предназначена для моделирования упругих и электромагнитных полей сейсмоэлектрического эффекта второго рода при возбуждении в одном из слоев горизонтально-слоистой среды механических колебаний точечным источником переменной силы.
1. Постановка задачи
Математическая модель сейсмоэлектрического эффекта второго рода в низкочастотном приближении может быть представлена в виде последовательности решаемых подзадач [1,2]:
- задача пороупругости;
- задача электрокинетики;
- электродинамическая задача.
Задачу пороупругости будем рассматривать в постановке Био [3, 4]. В соответствии с ней компоненты т тензора напряжений имеют вид
т. = 2^еи +Ъи (Х се -aM c)),
(1)
где е = divu ; q = -divw ; w = f (uf -u) ;
1
е.. = —
lJ 2
du. du. — + ■ J
dx. dx ..
компоненты тензора де-
формации; / - пористость; и - вектор смещения твердой фазы; и/ - вектор смещения жидкой фазы; Ъу - символ Кронекера.
Давление флюида определяется выражением
Р = -аМе + Мс. (2)
В выражениях (1), (2) постоянные а, М, Xс определяются соотношениями:
1 —
Хс = А + 20 + Я, аМ = — (0 + Я), м = —,
/ /2
О-/)
А -■
( ^ ^ 1 -f-tbL
K-sJ
V
K,+f
К-К» к,
-ILI'
Q =
f K л 1-f-±m_
K У
к
, R = fK
к к К = 1 - f--m + f_L.,
к J Kf
где ц - модуль сдвига скелета; Km, Ks, Kf -
модули всестороннего сжатия скелета, твердой фазы и жидкости соответственно. Уравнения движения имеют вид
дт,7 ч
1 + — =-ш2(рЛ +р ч),
дх.
дР
---+ Г/ = -ш2р - /ш У(ш),
дх,. 7
где рв = /р/ +(1 -/)р,; р/, р, - плотности жидкости и твердой фазы соответственно; У (ш)
- оператор Био; —, - объемные внешние силы, действующие на скелет и флюид соответственно. Для оператора Био будем использовать следующее выражение:
У( ш)=
1-1
. 0) М„
\ 1/2
2 ,
-1-
сэ
/ Л
где шь =---критическая частота Био;
а» к0 р/
Мь = 1-2 (зависит от формы пор); ада= 1-8 - извилистость пор, п - динамическая вязкость жидкости, к0 - проницаемость среды в стационарном поле (ю ^ 0).
После нахождения вектора из задачи пороупругости решается электрокинетическая задача. Пренебрегая обратным воздействия электрического поля на механическое движение (низкочастотное приближение), плотность стороннего тока находится из соотношения
,
(3)
где у(ш) = -¡шЬ(ш)У(ш) ;
Г - „ V1'2
Ь((й)=Ь0
1 . ю 2
1- г — -
ч % Мъ J £ = $ , 8/ - диэлектрическая проницае-
а„ Л
мость флюида; С - дзета-потенциал.
На последнем этапе решается электродинамическая задача для уравнений Максвелла
Ю Н = оЕ + )х го1 Е = / шд0 Н,
(4)
(5)
где Е и Н - напряженности электрического и магнитного поля; О = О-/ш8 ; О- проводимость; 8 - диэлектрическая проницаемость. Среды предполагаются немагнитными: ц„ = 4п ■ 10-7 Гн/м.
2. Расчет упругого поля сейсмоэлектрического эффекта
Пусть в круговой цилиндрической системе координат плоскости 2 = Н1 = 0 , 2 = Н., (Н] > 0, ] = 2,3,...,п, Нр > Н1 при р >1 ) разделяют среду Био на п +1 статистически однородных слоя с параметрами ак , Мк, X ^к , р Гк, рвк , Ук (ш), Уk(ш), Ок , причем
ск '
к=
2 < Н„
Н1 < 2 < Н2, Н2 < 2 <Н3,
п + 1, 2 > Н_
и будем называть слои при к = 1 - верхним полупространством (слой 1), к = 2 - слоем 2, к = 3 - слоем 3 и т. д., а слой при к = п+ 1 - нижним полупространством. Введем также обозначения = Н. - Н-, у = 2,3,..., п для мощностей 2,3,...,п слоев (к1 =», кп+1 =»). Предположим, что на одной из поверхностей 2 = Нт, т = 1,2,..., п расположен гармонический источник радиальных механических колебаний
Г' = ер — 5(р- а)8( 2 - Нт),
где — имеет размерность Н/м. Будем считать, что Г7 = 0 .
При таком возбуждении упругих колебаний задача о сейсмоэлектрическом поле осесимме-трическая, и в каждом слое искомые величины не зависят от координаты ф, а отличными от нуля являются компоненты и р, иг, чр, Ч2.
На влагопроницаемых поверхностях раздела 2 = Н., у = 1,2,..., п справедливы условия
сопряжения
и
и
w
р( 1
т( 1+1) | -т. ^ и=Н,+0 ^ \х=Н,
I = и(1) I
1г=Н,+0 \2=Н,-0 '
I = и (1) I
1г=Ну+0 мр 1г=Ну-0 '
I = Ч 1) I
1г=Н,+0 "г 1г=Н1 -0 '
) I = Т( 1) I
\г=Ну+0 \г =Ну-0 '
) I = Р( 1) I
\г=Н1+0 г \2=Ну-0 ' Г( 1)
-0 = -§ §(р- а).
Кроме того, координаты векторов смещения твердой и жидкой фаз, а также координаты тензора напряжений и давление исчезают на бесконечности.
Выберем в качестве искомых величин иг, и , Р, т22, трг и рассмотрим для них преобразования Ханкеля:
/•ад
и р=|ой ^тё£,
/•ад
V = 1 V ^^ £,
/•ад
и г =10 й£ ,
ад
w г =10 ^г£ , /•ад_
Р = |орш£Р¥£,
/•ад
Тгг = 1 ^о(£р)ё£ •
Вводя вектор-столбец
Х( к) =
Д к)
„( к)
„( к)
„( к)
„( к)
„( к)
1 2 (к) = -(к) Г к) = -(к)
4 ''"б .(к) = ■=( к)
где = и \ Х2 = ^ , ХЭ = Трг ,
х( к) = -(к) х( к) = -(к) X к) = р( к) Х4 _ 1гг , Х5 _ ир , Х6 _ Г ,
в каждом слое к получаем систему обыкновенных дифференциальных уравнений первого порядка с постоянными коэффициентам
ах
(к)
. = в( к) х( к
(6)
В( к)
=
третьего порядка;
0 Ь (к) 2
Ь( к) 0
нулевая матрица
Ь
(к) = 1
!,( к)
-ь -ъ
(к)
Ъ(к) Ъ(к) 2 иЭ
о
(к) д(к)
Ъ к)
о
,(к) _
Ч( к)
,(к) г,(к)
-ъ
-ъ
(к) г,(к)
(к) Ъ(к) ^11
,(к)
,(к)
ЪЭк) =-£, Ъ™=\1»к,Ъ™=т¥к,
Кк>=1/Ак,Ь^ = %(2цк-Ак)/Ак, ЬУ= а,/А,,
V ^
у
10
СОУ,
■ 3 2
1Ш3рл
,
Заметим, что
А 4 ^
Я* Рк
где V - скорость продольных упругих волн в скелете.
Матрица Био В( к) имеет шесть собственных значений х[к)=^к\ Цк)=\?, = -Х™, , , у которых
действительные части А,® = = ^у ^,
у^к) _ ^уЮ положительны и величины у^к] (к)
определяются выражениями
,
г-к
,
где
®2Р вк , Н»УА
ак =--1--
А, Ак¥к
К =
г'юКа
к к
-2-
соР/а
/шК,
■ 3 2
1(0 Рл
АД
ск =
М„
*к=.
Яе ёк > 0.
Таким образом, для решения задачи пороу-пругости требуется найти решение системы (6): координаты вектора Х(к) (г) непрерывны на всех влагопроницаемых поверхностях раздела, кроме поверхности г = Нт, на которой
хЭт+1) (Нт + 0) = х3т) (Нт - 0) - Га11 (£а) .
Кроме того, Х(1)(г) ^ 0 при г ^-ад и Х(и+1)(г) ^ 0 при г ^ад.
Предположим теперь, что на поверхности раздела г = Нк-1 + 0 известно искомое Х( к) системы уравнений (6), т. е. задано Х( к) (Нк-1 + 0).
Тогда решение Х(к) задачи Коши в слое к можно записать в виде
(к ),„д = в В( к)(г-Нк-1^( к)
Х( к)(г) = е
Х( к)(Нк-! +0),
В( к)( 2 -Нк-|) ,
где е к 1 - передаточная матрица 6-го по-
рядка.
Матрица В(к) имеет различные собственные значения , и по теореме Сильвестра
¿=1
гпе Па) _
' ПР'-Ч") №
Если известно поле Х(1) (0), то с помощью
В( к)( 2—Н к 1)
передаточных матриц е к 1 можно найти
поле в любом слое среды, в частности в слое (п + 1) решение уравнения (6) можно записать в виде
Х('!+1)(г) = е
_ в
,
где
В( 1\
Г =
вр =Пе" 1 , р < Я ;
0 -—а/Да) 0 0
Продолжая поле X с поверхности 2 = 0 в нижнее и верхнее полупространства и учитывая убывание поля на бесконечности, получаем систему алгебраических уравнений относительно компонент поля на поверхности 2 = 0 : 1. При 2 < т < п -1 -
W(1) Х(1)(0) = 0,
W(п+1)Вп Х(1)(0) = - W(п+:)В пт+1¥,
где
W(1) = w(п+1) =
с (1) с(1) с 2,1 4,1 °6,1
(1)
С( п+1) с °1,1 с3,1
( п+1) п+1) с ( п+1)
(к)
у - 1 -я строка матрицы В( ).
2. При т = п -
W(1) Х(1)(0) = 0,
w (п+1)В 2 х(1)(0) = - w (п+1)г.
3. При т =1 -
W(1) Х(1)(0-0) = 0, W (п+1)В 2 Х(1) (0-0) = - W(п+1)В 2 Г. (7)
Таким образом, в матричном методе решения задачи Био вначале определяются Х( к) (Нк), к = 1,2,..., п +1 (значения векторов Х(к )(2) на поверхностях раздела Н), кроме поверхности Нт, на которой находится Х(т+1) (Нт + 0). Если т = 1, то находится еще Х(1) (0 - 0). После чего Х( к) (2) определяются по формулам
Х(1)( 2) = еВ(1)2 Х(1)(0),
Х(1)(2) = еВ(°( 2-Н-1) Х(1- 1)(Н,_1),
при I = 2,3,.,п +1, за исключением значения I = т + 1, для которого
Х( т+1)( 2) = еВ(т+1)( 2-Нт) Х(т+1)( Нт + 0), I = 2,3,., п +1,
а если т = 1, то
Х(1)( 2) = еВ(1)2 Х(1)(0-0). (8)
3. Расчет электромагнитного поля сейсмоэлектрического эффекта
Электрокинетический ток определяется в
каждом слое к, к = 1,2,.,п +1 горизонтально-
слоистой среды по формуле = ук (о^уг^, а
„,(к) = „.(к) „(к) Т поскольку ^ = „ р „ 2 и относительная
скорость флюида не зависит от азимутальной координаты ф, то и возбуждаемое электрокинетическим током электромагнитное поле не зависит
1ех( к) р
1
ех ( к )
.. -. ■ • сх
от этой координаты. Ток = играет роль стороннего тока и в силу азимутальной симметрии возбуждает в пространстве компоненту Н^ напряженности магнитного
напряженности
поля и компоненты
Е(к) Е(к)
электрического поля.
Представим компоненты электромагнитного поля и электрокинетического тока в виде преобразований Ханкеля:
Н ()=| С Хиш $,
0
Е(рк) =| ерк) &/Др>/ $,
0
Е(к) =|е[к^рУ$ ,
/;(к)=| ]рк ^дрУ $,
гх( к)=| ]2к V0($р)^$.
Тогда из (4), (5) получаем систему уравнений для вектора-столбца
Z(k) = ||z(k) Zflf = fek) e(k)|f,
1 2 ф P
k = 1,2, ..., n +1:
dZ( k)
dz
■ = A (k) Z (k) + F (k),
(9)
7(k) _ _ Ak) (
где F( k )(z) = F(k) F(k ЧГ F{k) =- 7
em ^ ' em 1 em 2 ' em 1
F(k ^ = — j k), а матрица коэффициентов A
(k) em
системы обыкновенных дифференциальных уравнений первого порядка имеет вид
0
A(k) =
em
S2
---+ /га ц0 0
ст,.
Собственными значениями матрицы Aem являются k) = g k, g2k) =-? k, где Re ^ > 0.
После определения вектора w(k) в соответствии с формулой (3) можно записать следующие выражения для элементов вектора-столбца
Fem )(z):
Fem) = I ef,
Z(n+1)
■ о , z2n+1)
0.
(14)
Найдем решение в слое к. Пусть задан вектор-столбец Z(*~1) (Н^). Тогда Z(к)(г) является решением системы (9) и подчиняется условию
г«(я,_1) = ^-1>(я,_1). (15)
Решением задачи (9), (15) будет
Z( k)(z) = e
_ Aem>(z-Hk-tw(k-1)
HJ + Z(k) (z), (16)
где
Aem)( z-H^ 1) _
= I
C,( k)(z-Ht- 1^(k)
G(k
1 1
G( k) = 2 Ck 2C k 1 , G( k) = 2 2 Ck 2C k 1
2 2
Z(k)(z) =
7(k) 7(k)
- решение уравнения (9),
Р1 Р 2
удовлетворяющее условию
Z(ík )(Ик-1) = 0.
Решение Zlpk) (г) может быть найдено, например, методом вариации произвольных постоянных и иметь вид
Zpk )(z) = eck(z-H-■)g(k) + e
< k (z-Ht-1) У (k)
- g- +
F( k) = I e1'k)(z-Hk - ')f ,(k), k = 2,., n,
em / 1 l 7 7
F(n+1) = I
em
n+[)(z-Hn )f (n+1)
где fl( k) j j = 1,2, l = 1,2, .,6; k = 1,2,..., n +1;
fK Ja
(k) _ Y к
г'соГ,
f (k) = Yk^ b (k).
J2l ~ - 2l ' CT,,
(10)
b() - элементы матриц b( ), определяемых равенствами
b(1) = D(1) X(1)(0-0),
(11)
Ь((к) = Б((к)Х(к)(Як_1 + 0) (для к = 2,3,..., п),
Ь(и+1) = Б(п+1) Х( п+1)(Яп + 0). (12)
Функции 2[ к-1,2\к-1 непрерывны на поверхностях раздела г = И,, I = 1,2,...,п :
)(И1) = ^)(И1+1), 1 2 \И,) = 14)(И,+1), а на бесконечности при г ^ +оо
1(1) ^ 0 , ^ 0. (13)
где g (k) =
у (k) 311
+ I e
l=1
g (k )|l &21
')(z-HM)g( k)
(17)
g 2k) =
j( k) 312
j( k) 322
g(k) = gik)
# г l=1,2, .,6;
1 6 r /"(i:) - rr oll — ^ Zj
2C
к 1=i
ж
1 6 r f(i) -1- гт /'(i) Zj ¡те
^12
2C
к M
С,
g(k) = - Ct g(k) g(k) = Cl g (k) &21 - &11 ' &22 - &12 '
-(i) _ /1/ kJ21 61/ _'
l
у W - f(k) f if K .
, k = 2,3
V4 -
•2 Лк)
sv =■
? • • •?
n.
При k = 1 выражение (17) примет вид
Z(;4z) = eCiz g((1) + e-C'z g 21) +I e^'<1)z glu ,(18)
k
l=1
1=1,3,5
k
k
1=1,3,5
где
1 Г /"(1)
ГО
,
&5 + ^
а при к = и + 1 -
Х "+1)(2) = г-Н»)«(к)
g(к ) + е ^и+1(г Нп (к) +
^"+''( г-Н„ )~( к )
+ 1
(19)
где
/=2,4,6
«г=-
Г Лп+1) _ ~ г(п+1)
2С л+1 1^4,6 -
2С
— £
л+1 '=2,4,6
}+с
'л+1
Решение (16) позволяет продолжить электромагнитное поле с границы любого слоя в любой слой. Выполним, например, продолжение поля с поверхности г = Н1.
В слое 1 -
Ъ(1) (г) = е ^ (0) + Ър (г), (20)
в слое 2 -
Ъ(2)( г) = еА»( гН Ъ(1)(0) + Ъ(Д г)
(2),
ъ (2)( н 2) = а 2г(1)(0) + ър2)( Н 2), в слое Э -
ъ(Э) (7) = еА»( г-Нг) [А2г(1) (0)+гр,2)] + ъ(р3) (г), ъ(Э)(нЭ) = а32г (1)(0) +
где введены матрицы
А = ГК
I=1
рШ =
р
x а . ър-1)(нн) + ъ(рщ), . > э,
ър )(Н.), 1 < э
Продолжая этот процесс, получаем следующие выражения в слое к = 2,3,..., и + 1:
Ъ(к)( г ) = еАет)(г-Нк-1) Ъ(к-1)(Нк-1) + Ъ(рк)( г ) = ) [ А2-1Ъ(1) (0) + ¥(рк-1) ] + Ъ(рк) (г),
= еАетч г-Нк-1)[ Ак-1ъ(1)(
Ъ( к)(Нк ) = А 2 Ъ(1)(0) + ). (21)
В соответствии с этой формулой на поверхности г = Н получаем
Ъ(п)(Ни) = А 2Ъ(1)(0) + Е(
а в слое п + 1
(п+1)/_\ _ „А<т+1)( г-Ни )7(и)
Ъ(п+1)( г) = е
Ъ(и)(Нп) + Ъри+1)( г) =
= е
лет+"( г-Нп )Г А ^ (1)
[А2ъ(1) (0) + Ери)] + ъри+1) (г). (22)
Отсюда следует, что по известным значениям ъ(1)(0) можно вычислить электромагнитное поле в любом слое.
Для нахождения Ъ(1)(0) запишем , применяя выражения (16), (18):
1=13.5
Для выполнения условия (14) необходимо, чтобы
о « ъ(1)(0)+§ 21) = о.
Это влечет за собой два условия:
22» +1^(0) + ^ = 0,
^ +^г21)(0)+£22) = 0.
2б1
г
В силу условия |г(1) = (1) одно из этих со-
22 СТ1 12
отношений можно отбросить (например, второе) и записать
12» + 1^(0) = -. (23)
Аналогично в слое п + 1, применяя (13), (19), (22), получаем
2С
или
-А" -2 2£
л+1
И+1 Лп
и+1 /
+
Г1 - Л
л+1 у
2
(24)
~ о Р1 ГРг #11 •
После решения системы уравнений (23), (24) относительно Ъ(1)(0) определяются Ъ(к)(Нк),
е
I=3
к = 2,3, ..., п, а затем по формулам (20), (21) - искомые Z(k)(z), к = 1,2, ..., п.
4. Аналитическое решение задачи о сейсмоэлектрическом эффекте для однородного пространства
При разработке на основе рассмотренных алгоритмов программы расчета упругих и электромагнитных полей сейсмоэлектрическо-го эффекта второго рода необходимо наличие набора тестовых задач для отладки работы ком -пьютерной программы. С этой целью может быть использовано аналитическое решение задачи о сейсмоэлектрическом эффекте второго рода при возбуждении в однородной среде механических колебаний точечным источником переменной силы.
Пусть все евклидово пространство Я3 заполнено средой с параметрами слоя 1. Тогда из решения системы (7) с учетом (8), (11), (12) можно получить решение задачи Био для пространства в аналитическом виде
г < 0,
Х(2)(*) =
г > 0.
Аналитическое решение для электромагнитного поля сейсмоэлектрического эффекта второго рода в однородном пространстве, заполненном средой Био, имеет вид
т
/=1,3,5
-¿ГО
si /
О)
о 1/
ст,
+
. V X,mz II /V (1) ~(1)||Т „
/=1,3,5
Z(2)(z) = e-zlZ £ /=13 5
¿(1) £ia(1) Si/ Si/
+
/=1,3,5
где
g/0) =
/ = 1,3,5;
f1(1), f1 - определяются равенствами (10). Выводы
Предложенная математическая модель сейсмоэлектрического эффекта второго рода в горизонтально-слоистой среде Био может быть использована при разработки новых методов ГИС для интерпретации результатов измерений упругих и электромагнитных полей сейсмоэлектрического эффекта при зондировании околосква-жинного пространства.
Работа выполнена при финансовой поддержке компании Шлюмберже.
Библиографический список
1. Светов Б. С., Губатенко В. П. Электромагнитное поле механо-электрического происхождения в пористых вла-гонасыщенных горных породах в 2 ч. Ч. 1. Постановка задачи // Физика Земли. 1999. № 10. С. 67-73.
2. Губатенко В. П., Светов Б. С., Московский И. Г. Электромагнитное поле механо-электрического происхождения в пористых влагонасыщенных горных породах : в 2 ч. Ч. 2. Расчеты в горизонтально-слоистых средах // Физика Земли. 2002. № 2. С. 34-50.
3. Biot M. A. Generalized theory of acoustic propagation in porous dissipative media // The J. of the Acoustical Society of America. 1962. Vol. 34, № 9. P. 1254-1264.
4. SchmittD. P., Bouchon M., Bonnet G. Full-wave synthetic acoustic logs in radially semiinfinite saturated porous media // Geophysics. 1988. Vol. 53, № 6. P. 807-823.
Образец для цитирования:
Губатенко В. П., Московский И. Г. Матричный метод решения задачи о сейсмоэлектрическом эффекте второго рода в геологической среде Био // Изв. Сарат. ун-та. Нов. сер. Сер. Науки о Земле. 2016. Т. 16, вып. 4. С. 241-247. Б01: 10.18500/1819-7663-2016-16-4-241-247.
