Моделирование сейсмического поля в акустическом приближении двухфазных, иерархически неоднородных сред

Процессы разработки нефтегазовых месторождений связаны с движением многофазных многокомпонентных сред, которые характеризуются неравновесными и нелинейными реологическими свойствами. Реальное поведение пластовых систем определяется сложностью реологии движущихся жидкостей и морфологического строения пористой среды, а также многообразием процессов взаимодействия между жидкостью и пористой средой [1]. Учет этих факторов необходим для содержательного описания процессов фильтрации за счет нелинейности, неравновесности и неоднородности, присущих реальным системам. При этом выявляются новые синергетические эффекты (потеря устойчивости с возникновением колебаний, образование упорядоченных структур). Это позволяет предложить новые методы контроля и управления сложными природными системами, которые настроены на учет этих явлений. Таким образом, пластовая система, из которой необходимо извлечь нефть, представляет собой сложную динамическую иерархическую систему.

При построении математической модели реального объекта необходимо в качестве априорной информации использовать данные активного и пассивного мониторинга, получаемые в ходе текущей эксплуатации объекта. Решение обратных задач имеет огромное значение для нефтяной промышленности, поскольку нефтяной пласт относится к числу природных систем, не поддающихся прямым измерениям и наблюдениям в целом. Исследования последних лет показали, что в эволюции динамических систем играют неустойчивости, природу которых изучает теория самоорганизации или синергетика. Информацию об их проявлений в нефтяном пласте при его отработке можно только получить, используя данные мониторинга, чувствительные к его иерархической структуре.

Настоящая работа посвящена выводу интегральных уравнений двумерной прямой задачи для сейсмического поля в динамическом варианте. В работе произведен совместному анализу интегральных уравнений двумерных задач для сейсмического поля в рамках модели локальной иерархической неоднородности с пористым включением и чисто упругой неоднородности иерархической структуры в приближении, когда параметр Ламэ μ=0, как во включении, так и во вмещающей его среде. В этом случае динамиче- ская задача сейсмики может рассматриваться независимо для случая распространения продольной и поперечной волны. В данной работе будет рассмотрен первый случай для предложенной модели. Полученные результаты могут быть использованы для выбора критериев комплексирования сейсмических методов исследования сложно построенных сред.

1. Задача о дифракции звука в двумерной пористой влагонасыщенной неоднородности n-слойной среды

Эту задачу будем решать, используя подход, изложенный в работах [2 – 4]. Массовые силы Ф будем считать потенциальными и сосредоточенными в первом слое n-слойной среды. Плоскость ХОY совпадает с верхней плоскостью 1-го слоя, z =0. Ось OZ направлена вертикально вниз. Образующие двумерной неоднородности в виде цилиндра произвольного сечения S ₀ направлены вдоль оси ОY . При μ =0 в каждом из слоев S _i выполняется первое уравнение из системы уравнений для прямой динамической задачи сейсмики [2], преобразованное к виду:

, 2        . ,

А ф . + ^.ф. = - ²ⁿ f ⁽ M);

CT;

u = дгад ф ; f ⁽ M) = —— Ф.; 2 ПЛ.

где i=1,…n , Ф i = Ф при i =1, при i ≠ 1 Ф i =0. Волновое число в i -ом слое равно согласно [2 – 4]:

2 2   2σi k1i=k1 =ω λi ;                                          (2)

где ω — круговая частота, σ _i ,λ _i — плотность и коэффициент Ламэ i -го слоя n -слойной среды. Пусть в J -м слое n -слойной упругой среды находится пористое влагонасыщенное включение (вода или нефть). Согласно [5], чтобы оставаться в равновесии жидкая фаза должна испытывать одинаковое во всех точках образуемого порами странства гидростатическое давление p 2 . Это давление должно действовать и на щающую среду. Обусловленная им ее деформация должна сводиться к изменению ема фазы V ₁ и пор V ₂ в одном и том же отношении:

про-вме-объ- мак-

AV₂     1

= P2. ^V2   ^K0 ^{2 ,}

K0 — истинный модуль сжимаемости фазы. В том же отношении меняется и весь роскопический объем среды. Для покоящейся в пористом включении жидкости внутри включения система (1), согласно [5] после несложных преобразований записывается в виде:

, 2.

^ Ф _i + ^k _Xi Ф _i = ^{- 2 n} ( fi ( M ) + ^);

CTiK

u = grad ^ ; f ⁽ M ) = —— Ф ; ^ = ⁽1 - x ^- ^) p_? ;

i        2n/i- i                KQ

V ( M ), M eS ₀

0, M 4 S ₀

K = λ — модуль всестороннего сжатия при μ =0, χ — пористость.

Введем обозначение: _ɶ

k ( M } = k \[ — волновое число в слоистой среде S i , i = 1,...,n

и

K ( M )=

k 1 ji при M ^£ Sg k ( M ) при M ^ Sg

k2 =ω2σji;

1 ji λ _ji

Индекс ji обозначает свойства среды, внутри неоднородности S0 . В общем случае в произвольном слое Si или внутри неоднородности S0 уравнение (4) с учетом (2), (5) и (6) будет иметь вид:

^AP z + k^= ^{- 2 n (} fi ⁽ M ) + ^ );

Граничные условия в среде без разрывов заключающиеся в непрерывности вертикальной составляющей вектора смещения и компонентов тензора напряжений, согласно [2], на границах раздела L i имеют вид:

Ni ^д Ф i + 1 = _{0      ;} f дф ) = f дф) _;

5z      Sz       ze L' ’ ( dn ) ji Vein ) ja ’ уi (to ^i +фi)

^у i + 1^{( to 2} ^ i + 1 ^{+ Ф}i + 1) =⁰l

z e Li^;

у ( to ф + Ф + a p 2)

= ^{у ( to 2} ф + Ф ) ^- ji ^L

- ja

^{у ( to 2} ф + Ф ) |z= o = 0;

индекс ji – обозначает значения σ ,ϕ,Φ на границе неоднородности с внутренней стороны, jа – с внешней стороны границы неоднородности, которая расположена в j-ом слое, L-граница раздела слоя, с индексом i – со стороны i-го слоя, с индексом i+1 – со сторо-K ны (i+1)-го слоя, a=1—x-K0 , согласно (4).

Условия затухания на бесконечности согласно [2] имеют вид:

r grad 9 i = O (1),

r⁽ — ik_№ ) = o⁽ⁱ⁾ . o r

Пусть:

ɶ

9 i = T i — Ф 0 ,

где i=1,...j,ji,...,n, $ - потенциал нормального сейсмического поля в слоистой среде в отсутствие неоднородности: ϕ ji = ϕi и

A 0 . 7 2 0 n V

A p i + k _1ip i = ^{- 2n} fi ⁽ M );

^d Pi а^

d z      d z

= ⁰ z e Li- ;

с i ( to 2 ^ 0 + Ф0

с i ₊ i⁽ to ² ^Р _{+ 1} + Ф/ ₊ 1) = 0 z e Li .

0    ∂ϕ ⁰

На контуре неоднородности ϕ и _∂n непрерывны. Алгоритм вычисления нормального поля при дифракции звука в слоистой среде для произвольного источника возбуждения изложен в работах [6–7].

ϕ_i – потенциал аномального сейсмического поля, который, как легко показать с учетом (2),(4),(7) и (11) удовлетворяет уравнению:

Д p + K ²( M 'Р =- ( K ²( M ) - 1 ²( M ))( p 0 + a P2);

и граничным условиям:

d pi   d(pi + 1 _ _n        f d p ) _(d p )   .

0          ;                                ;

^dz      d z           i \ d n J ji (dn ) ja

2-σ _i ω ϕ _i

⁷ i + 1 ® Pi + 1

^{— 0} z e Li ^;

На контуре неоднородности:

[ ^^^ ] ja - [ ^^^ ] ji = ^(o" ja ^{- a} ji ⁾⁽ ^ j ° ^{+ a} P2)                         (15)

Функция источника сейсмического поля GSp (M,M ) определяется как решение следу- ющей краевой задачи [3–4]:

AGq . + k G^ .(M, M ) = -2пЗ(M - M ); Spi Spi и граничным условиям:

^{d G}Sp,i    ^{d G}Sp ,i + 1

dz        dz

= ⁰ zeLi ^;

CTto²Gjjp |z= 0 = 0;

2^„

⁶⁷ i ® ^GSp,i

⁷ i+1 ® ^GSp,i + 1 = ⁰ zeLi •

На контуре неоднородности

GSp ( M , M ⁰) и

d Gsp ( M , M ⁰) ∂n

непрерывны. Применим

формулу Грина [2] для функций ϕ _i ды при i^j :

и G Sp ⁽ M ’ M ) для каждого слоя n -слойной сре-

J («Pi

Li

dGSp , i дп

дф;

Go .-d- ) dl. Sp,i дп    i

_f^dGSp , i +1

J ⁽Фi +1---f---- д п

—G

^d ^i +1

^SP , i +1 дп

) dli +1

2 nq ) i ( M ⁰) при M ⁰ e Si 0 приМ ⁰^Sj

При i=j :

I ( Ф j ' "j

LJi       '"

d(pj

G^   ^j )dl •

SP,J dn    j

I ( 0 j + 1 Lj + 1

^{d G}Sp, j + 1 dn

5^7+1 „

G^    —J—) dl i +

Sp, J + 1    dn     J + 1

( - dGSp, j + I (Ф ja   PJ on c ja

^d ^ ja..

Go • —z—) dc • „ = Sp, j dn     ja

2 лф j ( M ⁰) приМ°eSj -S ₀ 0 при M ⁰ ^Sj - S o

C ja – контур неоднородности С (внешняя сторона).

Домножим выражения (18) и (19) на σ i и сложим их с учетом граничных условий, в результате чего получим:

°!  a®^

2n pVJa 0n

^{В ф}ja^r J GSP, j Bn ) dc "I

о ( M0)ф ( M 0) приМ0eSj -S0

приМ 0£Sj -So

T               M^\

Применим формулу Грина для функций ^ ji и GSp , j ⁽ M , M ) для внутренности области S 0 , с учетом (4), (5), (11), (12), (13), (14), (15), (16) и (17) получим:

k2" -kГ                   о ji 1 j fl P(M)GSp,j(M, M0)dTM

²ⁿ   S o

1 ~ ^dGSp, j ^J ( p ji ^J 2 п c ^J dn

∂ϕ _ji

— G_c  . —д— ) d,c =

Sp , j d n ’

( p ( M ⁰) + a p2) приМ ^O eS o 0 npuM ^O ^S o

Домножим выражение (21) на σ ji и сложим полученный результат с выражением (20). С учетом граничных условий на контуре неоднородности получим:

^ji ⁽ kxn kу)                    о         ^(Gia-^j0_£ 0     dGsp /

-----Tj--- L и ф ( M)G_Sp(j, M⁰)dT  + ■^{j   j} If Ф⁰(M)-^p^jdc

2n     _c                    M       2n              d n

SO                                         с

-

⁽c. ja " ”.ji\       6 _Ф

If GSp, j ^^{d -2n r J dn}

c

Воспользуемся равенством [2]:

( ^ ja — ^ji ) 2 n

_d 0                                    (22)

f GSp, j     ^dc= °' ⁽M^{0 )(}.Ф ⁽M^{0 ) +a} P2);

∂ n c

^°j^ »( , °( M )      -G_Sp , j⁸-4^M ) dc = J

^{2 n} c          dn            dn

⁽ ° ja —°" ji ) Ф °^{( M} °) при ^M ° ^S°

О                  (23)

О при M gS°

Тогда выражение (22) с учетом (23) можно переписать в виде:

( k ..—k                                          _CT.

j j II Ф ( M)GSp,MM , M ⁰) dT_M + -j- Ф ⁰( M ⁰)

-

^{2 n}     S₀        ^{p ,} j          M    ^о ji

^{( о} j- ^{- о} ji )

^о ji 2 ^П

4 GSp,j dp ^dc = ( Ф ( M ^{0) + a} P2) ^пРи^{м 0}^ ^S0 cⁿ

° ji ^(ki ц k 0                                0   0

-----j---j- ff ф(M)Gspj(M, MO)dT + ф0(MO) - o(M0)2д SC P J         M

-

^{( o} ja °ji ) , дф      /,,0.      „0   „

----о--- 4 GSp, i^rdc = Ф ⁽M ) ⁿP^uM ^ So

° ( M ⁰)2n _c       ^{d n                      0}

(24’)

Таким образом, решив интегро-дифференциальное уравнение (24) и определив распределение потенциала вектора смещений внутри неоднородности, мы можем, используя второе интегро-дифференциальное представление (24’), определить потенциал вектора смещений в любом слое, а затем, используя соотношение (1), вычислить распределение вектора смещений в любом слое.

• 2.    Сопоставление алгоритмов моделирования сейсмического поля для случаев упругого и пористого включения

Сравним полученные выражения с решением задачи дифракции сейсмического поля в рамках той же геометрической модели но с различными физическими свойствами включения. В работе [3] выписана система интегро-дифференциальных уравнений для случая упругого включения в n -слойной среде:

(kl ,7-kl                                      °й

• ^j, ^j If Ф( M)Gs_P(MM , M °) dT_M +     ф°( M °) -

• 2n    s°        p’j           M ° ji

-

^{( o} ja—° ji ) ° ji ²ⁿ

^1 GSp,j ^ ^dc =Ф (M ⁰⁾ при ^м °eS0 c ⁿ

°" ji ⁽ kx    ^kl /)                          о 0    0

• ^{1    j}' ^{1 j} II Ф ( M)GSp,j ( M,M°)dTM +Ф °( M °) - o ( M °)2n SC

-

^{( ct} ja ^ji ) , дф °^х     °f^ „ о

----0---^F GSp, j ^^dc = Ф ⁽M ) ^п Ри ^м ^ So o-( M °)2п c      ^{d n}                       °

Обозначения в (25) те же, что и в (24) и (24’). Сравнивая (25) и (24), легко заметить, что имеет место различие в структуре свободного члена в интегро-дифференциальном уравнении для внутренней задачи. Это, безусловно, скажется на различии в решении внешней задачи для этих двух моделей. Однако наличие пористого включения не приводит к изменению волнового числа для рассматриваемой задачи распространения продольной волны в рамках двух разных моделей среды, что может свидетельствовать о

неинформативности кинематических характеристик продольных волн для идентификации пористого влагонасыщенного включения.

• 3.    Моделирование дифракции звука

на упругой неоднородности иерархического типа с пористым включением

Рассмотрим задачу дифракции звука на двумерной упругой неоднородности с иерархической структурой, расположенной в j -ом слое n -слойной среды [4]. Если при переходе на следующий иерархический уровень ось двухмерности не меняется, а меняются только геометрии сечений вложенных структур, то аналогично (1) можно выписать итерационный процесс моделирования сейсмического поля (случай формирования только продольной волны).

Идею, изложенную в работе [4] для решения прямой задачи для двумерного случая распространения продольной волны через локальную упругую неоднородность с иерархической структурой, расположенной в J -ом слое n -слойной среды, распространим на случай, когда на L -ом иерархическом уровне окажется пористое влагонасыщенное включение.

σ

+ ^ja

<7

jil

P i-i (M 0 ) -

' k ₂_ k j jj P i ( M)G spj (M , M 0 )dT m

S _Cl

(j^ ^ Gspj а 21л jil             c

^c = P l (M ⁰ ), M ⁰ ∂ n

∈ S Cl

a _Pi ^{( k}J ^- kj) a ( M ⁰ )2n

jj P i ( M W^M , M ⁰ )d T M

+ P i⁰-1 ( M ⁰ ) -

^S Cl

^{( a} ja — ^а Jil )

a ( M ⁰ )2n ^f G ^Spj

^d l^l- dc = P | (M ⁰ ), M ⁰ ∂ n

∉ S Cl

где G sp i ( M , M ⁰ ) — функция источника сейсмического поля, она совпадает с функцией выражения (16), (17), ki jn = to ² ( ^ j^ 1 2 j_U ) - волновое число для продольной волны. В приведенном выражении индекс ji обозначает принадлежность свойств внутри неоднородности, ja – вне неоднородности, l =1… L-1 – номер иерархического уровня, гц = gradф l , ф l - потенциал нормального сейсмического поля в слоистой среде в отсут- 0         00

ствие неоднородности предыдущего ранга, если l=2…L, ϕl = ϕl-1 , если l=1, ϕl = ϕ , что совпадает с соответствующим выражением из [4]. Если при переходе на следующий иерархический уровень ось двухмерности не меняется, а меняются только геометрии сечений вложенных структур, то аналогично [4] можно выписать итерационный процесс моделирования сейсмического поля (случай формирования только продольной волны). Итерационный процесс относится к моделированию вектора смещений при переходе с предыдущего иерархического уровня на последующий уровень. Внутри каждого иерархического уровня интегро-дифференциальное уравнение и интегро-дифференциальное представление выписаны для потенциала, через которого выражается вектор смещений (26). Если на некотором иерархическом уровне структура локальной неоднородности распадается на несколько неоднородностей, то двойной и контурные интегралы в выражениях (26) берутся по всем неоднородностям. Если l=L, то внутри неоднородностей предыдущего иерархического уровня оказывается пористая влагонасыщенная неоднородность. В этом случае система (26) с учетом (24) и (24’) переписывается в виде:

^{( k}Qil -

2 n

k ' fl № ( M )G*,/ M , M 0 )d T M

S _{0 l}

+ ^ja M . (M 0 ) — σ jil

( ⁷ ja ."i ) fl G , ^dc = ( № ( M ⁰ ) + ap ₂ ), M ⁰ e S

<7 ,2n   J       on jil              c

"jil^ »^ J"I ^ ' ( M ) G spj (M ’ M 0 ) d T M + ^° 1 ( M 0 ) - ^ ^{( M )2 n   J} s₀,

• - (°" ja ^- " j.1 ) I G    SP1_ dc = ^ ( M 0 ) M 0 ^ S 1 = L

"( M ⁰ )2n      ^Sp,j d n        ^{l                01}

Если l = L +1 и на следующем уровне неоднородность снова упругая, то для дальнейшего продолжения итерационного процесса мы снова можем использовать выражения (26).

Заключение

В работе рассмотрена задача моделирования сейсмического поля (случай распространения продольной волны с учетом только сжатия и растяжения среды) в слоистой среде с включениями иерархической структуры. Построены алгоритмы моделирования в сейсмическом случае для 2- D упругой и пористой неоднородности. Исследован вопрос отражения физических свойств однофазности и двухфазности в системе уравнений решения прямой динамической сейсмики в частотном варианте. Представляет интерес с использованием полученных алгоритмов исследовать вопрос об изучении связи между тензорами напряжения и деформации на каждом иерархическом уровне и о возможном отклонении ее от обобщенного закона Гука. С другой стороны с увеличением степени иерархичности среды увеличивается степень пространственной нелинейности распределения составляющих вектора смещений, что предполагает исключение методов линеаризации задачи при создании методов интерпретации данных динамической сейсмики. Кроме того, усложняется процесс комплексирования методов, использующих волновые поля для изучения отклика среды с иерархической многофазной структурой. Эта проблема неразрывно связана с формулировкой и решением обратной задачи для распространения электромагнитного и сейсмического полей в таких сложных средах, что является нашей ближайшей задачей.