Моделирование кислотного гидроразрыва пласта

Во введении обоснована актуальность выбранной темы исследования, представлены цели и задачи работы, выделены положения представляющие научную новизну, сформулированы основные результаты работы и положе- ния, выносимые на защиту. Также выделена теоретическая и практическая значимость работы, приведены сведения об апробации результатов и основ- ных научных публикациях. Кратко изложено содержание диссертации по главам.
В первой главе представлен обзор технологий кислотного воздействия (КВ) и кислотного гидроразрыва пласта (КГРП), а также математических моделей и методов моделирования КВ и КГРП, известных в научной лите- ратуре.
Основными методами КВ являются: кислотная обработка, кислотная ван- на и КГРП. Кислотная обработка – наиболее распространнёный тип КВ. В карбонатных коллекторах она применяется с целью интенсификации добычи и увеличения приёмистости нагнетательных скважин. При КГРП в скважи- ну нагнетают кислотный состав при давлении достаточном для образования трещины разрыва. После снятия давления трещина не закрепляется и смы- кается, стимуляция обеспечивается неравномерным растворением породы, зачастую локализованном вблизи стенок трещины.
Большой вклад в развитие математических моделей и методов КВ внесли
Б. Г. Логинов, R.S. Schechter, R. Lenormand, H.S. Fogler, M. Quintard, Г.Т.
Булгакова, Р. Д. Каневская, К. М. Фёдоров; в развитие методов моделирова-
ния КГРП – A. Settari, D.E. Nierode, A.D. Hill, D. Zhu, В. Г. Салимов. В данной работе рассматриваются следующие реакции растворения
2HCl + CaCO3 → CaCl2 + CO2 + H2O, (1) 4HCl + CaMg(CO3)2 → CaCl2 + MgCl2 + 2CO2 + 2H2O. (2)
Кинетика реакций растворения карбонатной породы может быть представ- 9
лена в следующей форме:
−rk=Ak(bw(Xa−Xa))exp −RT , (3)
где rk – скорость поглощения единицы вещества кислоты на единицу поверх-
ности реакции, Ak – константа реакции, n – порядок реакции, bw – молярная
плотность водной фазы, Ek – энергия активации, Xa,Xeq – молярная кон- a
центрация и равновесная молярная концентрация кислоты в водной фазе. Параметры кинетики реакции, используемые в расчётах приведены в Таб. 1.
Таблица 1 — Параметры кинетики реакции используемые в расчётах.
eqn Ek
Минерал n кальцит 0.63
доломит 0.3 − 0.6
ððððHCl Ek
Ak ð2·ð·(ðððð) ð3
7.291 · 107 9.4 · 1011−3n
R (K) 7.55 · 103
11.32 · 103
Неоднородность пласта и многомасштабная структура фронта распро- странения кислоты в коллекторе ограничивают возможности моделирования КВ и предсказания результирующей продуктивности скважины. Режимы за- качки кислоты, оптимальные с точки зрения результирующей продуктивно- сти скважины, приводят к формированию узких, кавернообразных каналов – червоточин, которые имеют ширину, многократно превосходящую размер пор, и могут достигать в длину нескольких метров. Подходы к моделирова- нию червоточин могут быть разделены на несколько групп: на основе крити- ческого размера пор (R.S. Schechter), на основе модели капиллярных трубок (A.D. Hill), сеточные и фрактальные модели порового пространства (H.S. Fogler, G. Daccord), на основе чисел Дамкёлера (C.N. Fredd) и неравновесные численные модели различных масштабов (Р.Д. Каневская, M. Quintard).
Задача моделирования КГРП также включает в себя рассмотрение транс- порта кислоты по трещине. Одномерные модели распространения кислоты по трещине рассматриваются в предположении о несжимаемости жидкости в трещине и утечках в пласт согласно формуле Картера. В данной работе пред- лагается на основе двух- и трёхмерных моделей распространения кислоты по трещине, явно вычислить концентрацию кислоты и величину утечек в пласт. В рамках данного подхода распределение скорости в трещине предполагает- ся функцией давления, известной из решения задачи о течении в канале с
пористыми стенками.
Определение результирующей проводимости трещины КГРП представ-
ляет особенную трудность. Большинство моделей КГРП основываются на корреляции Нироуда-Крука, полученной на основе лабораторных экспери- ментов с разрезами кернов, которая определяет зависимость проводимости трещины от сжимающих напряжений пласта, его прочностных характери- стик. В подобных корреляциях также могут быть учтены оценки объёма и формы растворенной породы, влияние неровностей краёв трещины на её про- водимость, минералогическая и петрофизическая неоднородность её краёв. В данной работе предлагается подход для оценки проводимости на основе математического моделирования.
Во второй главе представлена математическая модель, используемая для моделирования КВ и КГРП. Фильтрация двухфазной водонефтяной мно- гокомпонентной смеси в неподвижной пористой среде подчиняется следую- щим законам баланса массы компонент:
∂ (msbwXi) + ∇ · (bwXiww) = qi, (4) ∂t
∂ (m(1−s)bo)+∇·(bowo)=0, (5) ∂t
∂ ((1 − m)bsk) = qsk, (6) ∂t
где s – водонасыщенность, 1 − s – нефтенасыщенность, молярные плотности фаз и скелета bw, bo, bsk считаются постоянными, скорости движения фаз ww, wo определяются согласно закону Дарси
wα =−Kkrα (∇p−ραg∇z), μα
(7)
где p – поровое давление, μα – динамическая вязкость фазы α, ρα – плотность фазы α, K – тензор проницаемости, g – ускорение свободного падения, krα – относительная фазовая проницаемость по фазе α, определяемая по формуле:
krα(m,sα)=m−m0k∗ +m∗−mkrα0, (8) m∗ − m0 rα m∗ − m0
где m0 – начальная пористость, krα0 – исходная относительная фазовая про- ницаемость фазы α, kr∗α = sα – предельно достижимая фазовая проница- емость, m∗ – максимальная пористость, достижимая при растворении. Эта
11

зависимость представляет собой средневзвешенное значение исходной фазо- вой проницаемости krα0 и предельно достижимой kr∗α. Связь пористости и проницаемости задаётся как
kl = eA(m−m0), (9) kl0
где l = {x,y,z}, k0,m0 – начальные проницаемость и пористость, коэффи- циент A – эмпирический параметр. Тензор проницаемости K предполагается диагональным K = diag{kx, ky, kz} и латерально-изотропным kx = ky, где kx,kz зависят от пористости в соответствии с (9). Молярные интенсивности qi,qsk определяютсяиз
q= r ν′′ −ν′ Ω (10) i kikik
k
где qi – молярная интенсивность образования (поглощения) вещества i-го
компонента, rk – скорость поглощения единицы вещества кислоты на еди-
ницу поверхности реакции, ν ′ , ν ′′ – стехиометрические коэффициенты i-го ik ik
реагента в уравнениях (1), (2) до и после k-ой реакции соответственно, Ω – площадь поверхности пор в единице объёма пористой среды
Ω m∗−m 2
Ω=sm∗−m , (11)
где Ω0 – начальное значение эффективной поверхности реакции при s = 1. Рассматривается течение двухфазной пятикомпонентной жидкости, вклю- чающей водную и нефтяную фазы α = w,o, а также нефть, воду, кислоту, соль и углекислый газ i = o,h2o,a,s,co2 в качестве компонентов. Водная фаза может как изначально присутствовать в пористой среде, так и зака- чиваться в виде раствора соляной кислоты. Также предполагается, что все продукты реакций (1), (2) мгновенно растворяются в водной фазе, нефтяная фаза включает лишь один компонент – нефть. Характерное количество угле- кислого газа, образованного в результате реакции незначительно и, обычно,
не оказывает существенного эффекта на свойства раствора.
Для простоты в дальнейшем также пренебрегается влиянием компонент- ного состава на плотности, вязкости, относительные фазовые проницаемости и скорость реакции (кроме концентрации кислоты в (3)), что, тем не менее,
может быть легко учтено в рамках модели. 12

В отличие от проппантного гидроразрыва пласта, при котором трещина закрепляется после снятия избыточного давления, при КГРП трещина не закрепляется и стимуляция достигается за счёт неравномерного растворения породы вблизи стенок исходной трещины. В данной работе под трещиной КГРП понимается именно стимулированная область с измененными вслед- ствие растворения свойствами шириной w вокруг исходной трещины гид- роразрыва. Мы не учитываем геомеханический фактор, рассматривая на- чальную трещину гидроразрыва фиксированной ширины w0 (Рис. 3а) на стадии закачки раствора. Поскольку после снятия избыточного давления незакрепленная трещина смыкается, при расчете остальных стадий КГРП принято, что w0 = 0.
Предполагается, что по трещине течёт несжимаемая жидкость, распреде- ление скоростей которой известно и задано аналитически на основе решения задачи Бермана о течении жидкости в канале с проницаемыми стенками. Уравнения баланса массы флюида и кислоты в трещине могут быть записаны в виде:
∇·v=0 (12) ∂Xa ∂ ∂Xa
∂t+v·∇Xa=∂yDef∂y , (13)
где Xa – молярная концентрация кислоты в трещине, Def – эффективный коэффициент молекулярной диффузии. Для случая течения по трещине нью- тоновской жидкости имеются следующие выражения для компонент скорости течения v = {vx, vy, vz} в плоскости трещины:
w 02 ∂ p 2 y 2
vx=−8μ∂x 1− w
, (14)
3 2y 1 2y 3
vy=vL2w−2w , (15)
w 02 ∂ p 2 y 2
vz =−8μ ∂z−ρg 1− w , (16) 0
где w0 – ширина исходной трещины, vL – скорость утечек жидкости сквозь пористые стенки.
Модели течения флюида в матрице и трещине сопрягаются на стенке
0

трещины в соответствии со следующими уравнениями:
(p) (f) (p) (p) ∂Xa (f) p| = p| , wwy| =vL, Xawwy| = XavL +Def
, (17) где верхние индексы p, f обозначают пористую среду и трещину соответствен-
но.
Для решения уравнений (4)-(6) была использована полностью неявная
схема метода конечных объёмов. Ниже представлены схемы в остаточной форме для уравнений (4)-(6) соответственно:
τ k n+1 Hn+1 = msXb |j,n+1 −τ qn+1 + n+1 Tn+1 Xb rα
(pn+1 −pn+1), (18)
ij i w j,n n+1 ij Vj
jβ i w μw
β jβ
j
β
τ k n+1
Hn+1 = m(1 − s)b |j,n+1 + n+1 Tn+1 b ro (pn+1 − pn+1),
oj oj,n Vj jβ oμo j β β jβ
Hn+1 = (1 − m)bsk|j,n+1 − τn+1qn+1, skj j,n skj
(19)
(20)
где индекс n+1 обозначает новый временной слой, индекс β – соседей ячейки j, Tjβ – проводимости между ячейками j и β, где проницаемость взята как среднее гармоническое с весами ширины ячейки в направлении ортогональ- ном их смежной грани, Vj – объём j-ой ячейки, τn+1 – шаг по времени. Ин- терполяция потокового члена в уравнениях (18), (19) осуществляется против потока.
В трещине для решения уравнений (12), (13) используется аналогичная полностью неявная схема метода конечных объёмов c двухточечной аппрок- симацией потоков. Схемы для уравнений (12), (13) могут быть записаны в следующем виде:
Hn+1 = τn+1 vn+1 ·njβ,
(21) (22)
– суммар- ная скорость потока между ячейками j и β, учитывающая конвективные и кондуктивные (диффузионные) механизмы переноса, njβ – вектор нормали к общей грани ячеек j и β, модуль которого равен площади этой грани. Для ап- проксимации межблочных давлений использовалась линейная интерполяция, для насыщенностей и концентраций – схема против потока. Модели трещины
jβ
Hn+1 = Xa|n+1 + τn+1 (Xav∗)n+1 ·njβ,
f1 j Vj
f2j n Vj jβ
β
где vjβ – конвективная скорость потока между ячейками j и β, v∗ jβ
β
∂y

и пористой среды сопрягаются на стенке трещины, полученные уравнения решаются в одной системе.
Представленные схемы линеаризуются и решаются итерационно методом Ньютона. Расчёт компонент Якобиана представляется трудоёмким, особенно в случае многофазной многокомпонентной фильтрации. В работе использует- ся функционал библиотеки ADOL-C, которая предоставляет удобный интер- фейс автоматического дифференцирования, включая расчёт разреженного Якобиана, основанного на алгоритмах раскраски графов.
На этапе решения производилось масштабирование уравнений в системе, исходя из размеров ячеек. Это позволило снизить обусловленность матрицы и уменьшить общее время расчёта за счёт возможности расчёта с большим временным шагом на начальном этапе. Параметры масштабирования были подобраны именно с точки зрения уменьшения обусловленности получаемой матрицы на первой итерации первого временного шага. Для решения полу- ченных систем линейных уравнений были использованы предобуславливате- ли ILU(0,1), ILUT совместно со стабилизированным методом бисопряжённых градиентов.
В третьей главе проведено моделирование соляно-кислотной обработки (СКО) в круговом однородном изотропном пласте с начальной проницае- мостью k0 мощностью h, рассчитана продуктивность после обработки по формуле
PI≡
q2πk0h drre
= , s=k0 −ln , (23)
re
∆p μlnre+s
rw rw
где re, rw – радиусы удаленного контура питания и скважины соответственно, s – скин-факторa.
В отличие от выполненных ранее работ, где рассматривались кинетиче- ская модель реакции первого порядка, в данной работе проведены расчёты СКО с кинетикой реакции дробного порядка в соответствии с результатами лабораторных исследований.
Соотношение между скоростью реакции и скоростью закачки характери- зует число Дамкёлера
Da = rkΩ0V0, (24) q
rk(r) r w
где rk – скорость реакции, определённая в (3), Ω0 – исходная площадь ре- акции, q – скорость закачки, V0 – характерный объём задачи (в расчётах считается равным объёму закачанного раствора).
Выполнен анализ чувствительности эффекта обработки к скорости закач- ки для случаев кальцита и доломита. Зависимости результирующего скин- фактора от скорости закачки имеют минимум, который соответствует ско- рости закачки, при которой обработка обеспечивает наибольшую продуктив- ность скважины. Этот результат подверждается многочисленными лабора- торными и численными экспериментами.
Однако в большинстве работ этот оптимальный режим растворения ас- социируют с распространением небольшого количества узких червоточин в неоднородной пористой среде. Тогда как в данном исследовании рассматри- вается крупномасштабная модель СКО однородного пласта, а эффект опре- деляется оптимальным соотношением между скоростью закачки и скоростью растворения.
Также было исследовано влияние кинетики реакции на продуктивность скважины после соляно-кислотной обработки коллекторов сложенных каль- цитом или доломитом. Глубина проникновения кислоты в пласт для этих двух случаев разная, что приводит к различиям в результирующих распре- делениях пористости и проницаемости. Это объясняется тем, что скорость реакции соляной кислоты с кальцитом на несколько порядков превосходит скорость реакции с доломитом. Также можно отметить, что профили про- ницаемости, соответствующие различным порядкам реакции (для кальцита), значительно отличаются. Более компактному профилю (n = 0.2) соответ- ствует более быстрая реакция за счёт более низкого порядка (основание < 1). В четвёртой главе исследовано кислотное воздействие в масштабах неод- нородного керна. Для учёта этой неоднородности было сгенерировано поле проницаемости в соответствии с логнормальным распределением с парамет- рами kG = exp ⟨ln k⟩ = 5 , σln k = 0.5, гауссовой вариограммой с корелляцион- ной длиной λ = 0.02 . Исходная пористость при этом оставалась постоянной m0 = 0.09. Была рассмотрена прямоугольная область размерами 0.4 м и 0.1 м толщиной 0.1 м. В изначально нефтенасыщенную область закачивался 15%- ый раствор соляной кислоты через нижнюю грань с постоянной скоростью. На противоположной грани поддерживалось постоянное давление, на боко- вых - условие отсутствия перетока. Производилась закачка фиксированного объёма кислоты Va = 10−3 м2, после чего закачка останавливалась, и про- должался расчёт стадии остановки, которая завершалась после отработки оставшейся кислоты. Рисунок 1 — Динамика распространения червоточины в керне при посто- янной скорости закачки в разные моменты времени: 1 – 0.028 час, 2 – 0.12 час, 3 – 0.3 час, 4 – 0.98 час. Представлены кон- центрация кислоты (a), концентрация соли (б), пористость (в), проницаемость и линии тока (г). В некотором диапазоне скоростей закачки растворение происходит с фор- мированием червоточин. На Рис. 1 продемонстрирована динамика распро- странения червоточин во времени при постоянной скорости закачки равной q = 4.17м2/час в разные моменты времени 1 – 0.028 час, 2 – 0.12 час, 3 – 0.3 час, 4 – 0.98 час. Представлены концентрация кислоты (a), концентрация соли (б), пористость (в), проницаемость и линии тока (г). В этом расчёте был рассмотрен керн, состоящий целиком из доломита. При данной скоро- сти закачки растворение происходит в режиме формирования доминантной червоточины. Расчёты проводились до момента прорыва керна. При фиксированном объёме раствора Va была исследована зависимость динамики растворения керна от скорости закачки. В рассматриваемом диапа- зоне скоростей закачки наиболее эффективное растворение (самый быстрый прорыв керна) происходит при наиболее высокой скорости закачки, что, в свою очередь, соответствует наименьшему числу Дамкёлера. Рис. 2а иллюстрирует рост длины червоточины с увеличением закачан- ного объёма кислоты для разных чисел Дамкёллера и скоростей реакции, характерных для доломита и кальцита. Показано, что в обоих случаях име- ется оптимальный темп нагнетания (число Дамкёллера), при котором осу- ществляется самый быстрый рост червоточин. В представленных примерах Daopt ≈ 75 для доломита, Daopt ≈ 26 для кальцита. Эти данные качественно отражают эффекты, полученные в экспериментах.