УДК 622.276
Н. А. Шевко
Пермский государственный технический университет
ГИДРОДИНАМИЧЕСКОЕ МОДЕЛИРОВАНИЕ ОКОЛОСКВАЖИННЫХ ЗОН ПЛАСТА
Разработана детальная численная модель околоскважинной зоны пласта, базирующаяся на несогласованных нерегулярных неортогональных разностных сетках различное топологии и позволяющая детально описывать сложные пространственные течения потоков флюидов с учетом геометрик стволов скважин и неоднородности првекважинных зон.
Гидродинамическое моделирование процесса фильтрации флюидов в околоскважинных зонах пласта (ОЗП) с детальным учетом их неоднородности повышает точность определения расчетных показателей работы добывающих скважин. Особенно это актуально при прогнозировании результатов различных геолого-технических воздействий (ГТВ) на пласт и ОЗП, проводимых с целью регулирования и повышения эффективности разработки залежей. Уровень достоверности при прогнозировании показателей эффективности ГТВ увеличивается при применении постоянно действующих моделей (ГЩМ) залежей за счет более полного учета имеющейся геолого-промысловой информации на стадии выбора и планирования технологий ГТВ. Однако использование моделей пластов в рамках существующих ПДМ, основанных на традиционных конечно-разностных методах, не позволяет детально учесть геометрию ствола скважины, неоднородность ОЗП, сложность околоскважинных потоков и особенности изменения гидродинамического состояния ОЗП в результате ГТВ.
Внутренний диаметр стенок скважин много меньше размеров ячеек, используемых в моделях пластовых систем, что создает большие трудности при гидродинамическом моделировании ОЗП. Это обусловило создание специальных разностных схем с моделированием скважины в узле (блоке) разностной сетки. При этом параметры потока на стенке скважины определяются через параметры потока в блоке, интегрально характеризующего некоторую область возле скважины. Приток к «точечкой» скважине по ребрам сетки (через грани блоков) при таком подходе неадекватно описывает приток к реальной скважине (не учитывается, например, радиальность притока флюидов). В области разработки этих схем можно выделить два основных подхода, учитывающих особенности течения возле «точечной» скважины: введение поправочных множителей в коэффициенты разностных уравнений, относящихся к скважине, и локальное сгущение сетки в окрестности скважины.
Первый подход, названный CCF-метод (Completion Connection Factor), развитый в работах Г.Г. Вахитова, А.А. Писмана, В.Г. Карея, A.M. Каземи, основан на нахождении связи между забойным давлением скважины и
давлением в ячейке численной модели и на предположении, что течение в области возле скважины одномерное радиальное стационарное однофазное. Однако такие предположения не всегда выполняются (например, предположение о стационарности и однофазности потока не выполняется при высоких скоростях изменения насыщенностей в околоскважинной области, а предположение о радиальном притоке - при наличии неоднородности околоскважинной области), что приводит к существенным ошибкам при определении расчетных значений обводненности, газового фактора и дебита скважины.
Альтернативным первому подходу является LGR-метод (Local Grid Refinement), основанный на предположении, что погрешность любой дискретной схемы прямо пропорциональна размеру использованных блоков (О. К. Розенберг, В. В. Освальдо и др.). В этом случае декартова сетка для залежи используется совместно с «произвольной» (радиальной, эллиптической) детальной ортогональной сеткой в околоскважинной области, которая расположена в одном или более блоках декартовой сетки. Детальные модели скважин, основанные на LGR-методе, используются многими исследователями, в том числе в зарубежных и отечественных коммерческих программных продуктах по гидродинамическому моделированию, в следующих случаях:
моделирование локальных эффектов (конусобразование воды и газа, обводнение);
моделирование горизонтальных скважин; получение скважинных псевдофункций;
получение корректного распределения давления возле скважины; детальный учет изменения каких-либо свойств в ОЗП;
проведение детализации исследуемой области без нежелательного измельчения глобальной сетки.
LGR-метод дает более реалистичное представление о геометрии течения и происходящих процессах вблизи скважины, однако при его использовании возникает ряд проблем:
необходимо использовать другую систему координат и метод дискретизации дифференциальных уравнений, отличный от метода, используемого для регулярной глобальной ортогональной сетки, что усложняет алгоритмы получения потоковых коэффициентов между узлами (блоками); нужно использовать простые способы сопряжения сеток различной геометрии и различных размеров, что затруднительно при сложных траекториях и форме стенки скважины;
наличие блоков малых размеров вблизи скважины существенно ограничивает временной шаг, необходимый для сходимости уравнений или обеспечения заданной точности решения в задачах многофазной фильтрации;
при большом числе скважин может быть получено слишком большое число блоков, что существенно увеличит размерность задачи и затруднит вычисления.
Для создания модели околоскважинной зоны пласта как объекта воздействия необходимо использовать подход, который позволяет более
эффективно описать геометрию ствола скважины, неоднородность ОЗП, особенности изменения гидродинамического состояния ОЗП в результате ГТВ. Это возможно при использовании разностных схем на нерегулярных разностных сетках (A.A. Самарский). На необходимость отказа от регулярных сеток и применения методов типа конечных или граничных элементов указывалось многими исследователями. Данный подход не имеет части проблем, связанных с простой геометрией блоков сетки в LGR-методе, а именно: для дискретизации дифференциальных уравнений (ДУ) на сетке с различной геометрией блоков используется один метод, учитывающий форму, размеры и неортогональность разностных ячеек; при построении гибридной сетки можно согласовывать ее геометрию с ожидаемой конфигурацией течения и формой стенки скважины и пласта, что позволяет меньшим числом блоков описать особенности фильтрации при заданной точности; сопряжение сеток различной геометрии и различных размеров производится значительно проще в силу большего произвола в выборе формы и размеров блоков локальной сетки.
Сложность решения задач с подобными сетками связана со сгущением их в окрестности скважины, а также с произвольной их геометрией: более сложный метод дискретизации ДУ (получения потоковых коэффициентов); наличие нерегулярной структуры разряженной матрицы при решении системы линейных алгебраических уравнений и применение специальной схемы ее хранения; применение высокоэффективных методов решения систем нелинейных ДУ многофазной фильтрации.
Эти недостатки ограничивали широкое применение метода, однако он применялся многими исследователями для моделирования комплексных (сложных) резурвуаров в двумерной двухфазной постановке (О.О. Чавент, A.B. Кохен, М.П. Джефри. Л.В. Юмарт, Д.Р. Форсайт, A.A. Дарлов и др.), в однофазной трехмерной (В.В. Гурчан), трехмерной многофазной постановке (A.C. Фунг, В.О. Хайберт, В.II Йотов). Реализация современных эффективных численных методов и использование высокопроизводительной вычислительной техники позволили применить данный подход для создания модели околоскважинной зоны пласта, рассматриваемой в настоящей работе.
Для математического описания большинства физических процессов, происходящих в нефтяных залежах при их разработке, применима модель трехмерной изотермической трехфазной трехкомпонентной нестационарной фильтрации сжимаемого флюида в неоднородной анизотропной по проницаемости деформируемой пористой среде, реализованная во многих коммерческих программных продуктах но гидродинамическому моделированию (ECLIPSE, MORE, VIP, FRAGOR, LAURA) - модель «нелетучей» нефти, полученная из обобщенных уравнений фильтрации Маске га-Мереса при пренебрежении растворимостью отдельных компонентов ,! фазах. В рамках модели «нелетучей» нефти возможен учет гравитационных, капиллярных и вязкостных сил, зависимости свойств фаз (плотность, вязкость, доля компонентов в фазах) и пористой среды (абсолютная и относительная проницаемость, пористость) от давления в фазах и насыщенностей фазами пористой среды. Данная модель фильтрации используется в модели ОЗП:
V • \у„ = т ° (т„р. 5»).....£>; В О х (О, Г],
«о = ~кКУ{ра
V• >у„ = т ° {т РЛ)-в ^ х
V• уу= т ° (т^лч+%рл)-ег-ец в ах (о,г].
Система дополняется следующими соотношениями:
--/>!> = Л (-О,
и уравнением состояния породы и фаз
К = К(.т.)■ к'(ро), ка - к(Л'а), и, =цО?„), та = т(ра), Р„ »р(р„), = ■'*(/>,,)■ Граничные условия на стенке добывающей Г ^ и нагнетательной Г. скважин:
Рп = Рг" на при на
К-п) = ?Г наГ^.Г^,
(\у; ■ п) — о',: на Г"„г'оа-, начальные условия ра = р™, 5„ = 5""', 5 = .
Детальное описание околоскважинных особенностей посредством разностных ячеек малого размера обусловливает чрезвычайно большие временные и вычислительные затраты, связанные со снижением устойчивости численных схем при наличии в области моделирования ячеек с малыми размерами. Обычное применение в этом случае неявных схем решения (1МР1Л51Т-методов), являющихся «безусловно» устойчивыми, сопряжено со значительными вычислительными затратами при обращении нерегулярной разрешающей матрицы с блочными элементами, с наличием ограничений на временные шаги из соображений точности аппроксимации по времени и обеспечения устойчивости итераций по нелинейности в связи со значительными нелинейностями коэффициентов уравнений при временных шагах, значительно превышающих условия устойчивости явных методов. Это существенно снижает эффективность неявных схем. Однако, учитывая тот факт, что изменение давления в ячейках сетки за временной шаг, ограниченный по условию устойчивости численной схемы, незначительно по сравнению с изменением насыщенностей в этих ячейках, имеет смысл разделить временные шаги для уравнений по давлению и насыщенности. Подобные подходы применялись на всем протяжении развития математических методов расчета процессов фильтрации при ограниченных вычислительных мощностях. Поэтому схема пространственно-временной дискретизации математических уравнений модели ОЗП основана на модификации явного (¡МРЕБ) и неявного последовательного (5ЕС>) методов (рассмотренных в работах Г.Б. Азиза,
Л.А. Сеттари, Б.В. Шалимова, М.М. Максимова, Л.II. Рыбицкой, В.А. Рождественского и др.), позволяющих разделить решение уравнений для давления (общего потока) и насыщенностей (долей фаз в патоке). Проведенные численные исследования подтвердили эффективность такой модификации численной схемы МРЕБ и БЕС^-методов для задач фильтрации в околоскважинном пространстве.
Для дискретизации трехмерной области вблизи скважины с некоторыми особенностями (трещина ГРП, перфорационные отверстия) автором разработаны специальные алгоритмы, позволяющие разбить область на тетраэдры, четырехугольные и треугольные гексаэдры. В общем случае, например при многозабойной скважине и сложной траектории стволов, на основе метода Делоне разработан способ разбиения ОЗП на тетраэдры (рис. 1). Реализация алгоритма дискретизации в этом случае сводится к следующему:
1) задается первый набор узлов, описывающих внешние границы области;
2) по данному набору узлов проводится первая (граничная) дискретизация;
3) задается окончательный набор узлов, учитывающий геометрию потоков флюидов и стенки ствола скважины;
4) по этим точкам дискретизации окончательно строятся тетраэдры. Такой подход позволяет строить экономичную дискретизацию со сгущающимися разностными сетками возле скважины.
Рис. 1. Сшивка сетки ОЗП скважины сложной траектории с регулярной сеткой пласта
Для построения разностных аналогов дифференциальных уравнений на неортогональных нерегулярных сетках существуют несколько обших методов, таких как вариационный (метод Ритца), вариационно-проекционный (метод Галеркина, метод наименьших квадратов, метод коллокаций и др.).
вариационно-операторные (метод опорных операторов, разработанный A.A. Самарским, A.B. Колдоба, Ю.А. Повещенко) и проекционно-сеточные методы, которые позволяют строить консервативные разностные схемы. Наиболее экономичным с точки зрения скорости и простоты вычисления элементов матрицы разностного аналога ДУ является метод опорных операторов. Поэтому консервативный разностный аналог ДУ для случая многофазной фильтрации строится (H.A. Шевко) на основе метода опорных операторов с использованием так называемого смешанного конечно-элементного метода (mixed finite element metod). В качестве опорного оператора выбран оператор GRAD, который аппроксимируется непосредственно. Для получения определяемого оператора DIV используют соотношения векторного анализа. Для рассматриваемой системы уравнений с учетом неортогональности и возможной несогласованности сеток (несоответствия узлов сеток на поверхности соприкосновения элементов из-за различной топологии или детальности) численный алгоритм (для lMPES-метода) будет иметь следующий вид. Пусть
(v,,v,,v,),vt -контравар. комп, непрерывные на я(ш,,.)|.
W, '
puAr, =а,ае R|,
е V,
(р,Г!' рС > Р^о Рг»,рг,,я,, е IVш.
1. По заданным на и-м временном слое рассчитывается общее давление на п Н-м слое с шагом, равным с1ТГ/. по уравнениям для давления:
(рГ1 , V ■ и)- {рг„ , и • И) - , и - п)г. - (Я;,, и • п)л , (\у", ,у|= V• и,¥бУл,.,
G ", е) = (z, V • е), eeVM,
+ (V • G", (У) + (qсо), Ф e Wllk,
(wu-n,£>) r,=(qu,0)rX, 0^Whk,
для всех ^еГ,пГг,
где
£2 л(т,, )еГ
2. Расчет значений насыщенностей на конец у-го шага размером Л5 п-го временного шага, используя давление р^1 для расчета общих потоков
уупт+1 = к
и насыщенностей для расчета долей фаз в потоке по уравнениям для водонасыщенности
(<;'>")= (/Г'-^Т1, и)- - - и), и газонасыщенности
{Сь/ , ф)= (V ■ , ф)-{/3^тРв,ф)~ Д гр,е, ф)~
Расчет насыщенностей производится до достижения равенства суммы шагов по насыщенности и шага по давлению
V
3. Проверка допустимости ошибки материального баланса компонентов с использованием давлений и насыщенностей =5'^' при соблюдении (1). В случае превышения допустимой величины ошибки шаг (б(/",/л+|) считается недопустимым, результаты расчета отбрасываются и величина нового шага уменьшается пропорционально величине превышения фактической ошибки над допустимой. При благоприятном исходе осуществляется переход на следующий временной шаг.
При аппроксимации ДУ при расчете потоков в неоднородной среде коэффициент т и компоненты К брались средневзвешенными по объему, фазовая проницаемость бралась «вверх по потоку».
Решение системы линейных алгебраических уравнений осуществляется градиентным методом (ЖТОМШ (усеченный метод ортонормированного
минимума) с предобусловливанием методами Irregular Nested Factorization или ILU.
Рассматриваемая модель реализована в программном продукте WELLSIM. Для проверки правильности программной реализации алгоритмов однофазной и многофазной фильтрации использовались известные аналитические решения стационарной и нестационарной фильтрации. Условные обозначения:
а - индекс, обозначающий нефтяную, водную и газовую фазы (о, w, g): K(.v,), т'(х,) - полный тензор абсолютной проницаемости и коэффициент
пористости породы в точке пласта с координатами Х!; к\р) - множитель к коэффициенту проницаемости как функция давления; ра, кп, р;, Liit. т а , S , /' - соответственно давление, относительная проницаемость, плотность, вязкость, относительная пористость, насыщенность и доля к-го компонента в а-й фазе; 1 ■ время; Т - продолжительность временного интервала; 0'с мощность источника (стока) к-го компонента, Q" > О - источник, QI <0 -сток; капиллярное давление между фазами а и В, используем р[ео, р1 ; g -
ускорение свободного падения; г - глубина залегания пласта; Q - исходная
область; Ol: - отдельная ячейка сетки, со..у - узел, ад...... грань, XR - ребро, где
индексы Е, N, В, R определяют глобальный порядок нумерации
соответствующих дескрипторов сетки; <о(ПЕ) ..... узлы, принадлежащие ячейке
Qg, ю„(Ор) - конкретный узел, где индекс п определяет локальный порядок-нумерации узлов и т.п. Д* ~{kaPj'a)l' Р* - подвижность компонента к в а-й фазе; f"'Lv = 4 )Аг (р.Г' ) - доля воды в потоке,
/Г!'У ~ - доля газа в потоке,
4 = Ы;: + сл°е + вц + в XI)", о" = (aqz + о°0 + bq* + bq* )",
где верхний индекс и - /2-й временной шаг; Г,, Г, - поверхности соприкосновения несогласованных (гибридных) сеток;
С« =т{т'рр*Х С« Со = "rWpP'o + т'^Ро)8"'
Д/ ¿S.t /\t
tf = ft («„А, г1. Cf = [п'р {pjs У + < (p'jf + рГГ/ к. ^ ^к^лГ с? =^пРр; + cs/
P.,, =- с"т/с!:, > fir = АС~ ссо - BC!i" - вс~х. PSo = cg°r/cga,
= с;<7 /с*,. A = i/c:> , в = i/csg<, С = fi - с* /с-;-1let'.
Рассматриваемая модель реализована в программном продукте WELLSIM, правильность программной реализации алгоритмов которой проверена аналитическими решениями однофазной и многофазной, стационарной и нестационарной фильтрации.
На основе созданной модели околоскважинной зоны пласта стало возможным решение следующих прикладных задач: приток флюидов к скважинам сложной траектории; приток к многозабойным скважинам и вторым стволам; приток к трещине гидроразрыва; приток к перфорационным каналам; приток при детальном учете неоднородности пласта и прискважинной области; проверка результатов интерпретации гидродинамических исследований скважин.
Получено 08.07.03 УДК 622.276 A.A. Щинанов Перм НИПИн ефть
МОДЕЛЬ ДВУХФАЗНОЙ ФИЛЬТРАЦИИ В ДЕФОРМИРУЕМОМ ТРЕЩИНОВАТО-ПОРИСТОМ ПЛАСТЕ
Представлена модель двойной пористости/проницаемости, которая позволяет моделировать высокоскоростную фильтрацию по системе трещин и учитывает динамическую деформацию трещинно-порового коллектора, возникающую при изменении пластового давления. Классическая модель единичной пористости сопоставлена с моделью двойной пористости/проницаемости. С помощью математического моделирования определены некоторые особенности фильтрации в трешиновато-пористой среде и влияние обмена между матрицей и трещинами на скорость фильтрации.
В работе рассматривается модель двухфазной фильтрации в среде с двумя видами пустотности, которая известна также как модель двойной пористости/проницаемости [1-3]. Базовой для данной модели является модель фильтрации в среде с одним видом пустотности (единичной пористости), широко известная как модель фильтрации в пористой среде [4].
Уравнения течения однородной жидкости в трещиновато-пористой среде с двумя видами пустотности были сформулированы Г.Н. Баренблаттом и др. [5] исходя из континуального подхода (условия непрерывности). По Баренблатту, обе среды - система трещин и блоки пористой матрицы - рассматриваются как две сплошные среды, вложенные одна в другую, причем параметры движения жидкости и среды определяются в каждой математической точке. Уравнения движения и сохранения массы записываются независимо для каждой среды. Переток жидкости из одной среды в другую учитывается введением функции источника-стока в уравнениях сохранения массы. Подход Г.Н. Баренблатта был распространен на случай многофазной фильтрации X. Каземи [2].