На правах рукописи

СОКОЛЮК Любовь Николаевна

СОВЕРШЕНСТВОВАНИЕ МЕТОДОВ

РЕМАСШТАБИРОВАНИЯ В ГИДРОДИНАМИЧЕСКОМ

МОДЕЛИРОВАНИИ ПЛАСТОВЫХ СИСТЕМ

Специальность 01.02.05 – Механика жидкости, газа и плазмы

АВТОРЕФЕРАТ

диссертации на соискание ученой степени кандидата физико-математических наук

Тюмень-2011

Работа выполнена в Тюменском государственном университете и в Тюменском филиале Института теоретической и прикладной механики им. С.А. Христиановича СО РАН

Научный руководитель: доктор физико-математических наук Родионов Сергей Павлович

Официальные оппоненты: доктор физико-математических наук, профессор Федоров Константин Михайлович кандидат физико-математических наук Богачев Кирилл Юрьевич

Ведущая организация: Московский физико-технический институт

Защита состоится «11» февраля 2011 г. в 1400 на заседании диссертационного совета ДМ 212.274.09 в Тюменском государственном университете по адресу: 625003, г. Тюмень, ул. Перекопская, 15а, ауд. 118.

С диссертацией можно ознакомиться в библиотеке Тюменского государственного университета.

Автореферат разослан « » декабря 2010 г.

Ученый секретарь диссертационного совета к.ф.-м.н., доцент Мусакаев Н. Г.

-2

ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ

Актуальность. В настоящее время при моделировании разработки нефтяных месторождений широко применяются трехмерные численные гидродинамические модели многофазной фильтрации. Гидродинамические модели создаются на основе геологических моделей нефтяных месторождений, которые в некоторых случаях могут содержать настолько большое количество ячеек разностной сетки, что время расчетов становится неприемлемым для практического использования модели. Поэтому становится необходимым уменьшить количество ячеек путем объединения мелких ячеек геологической сетки в крупные ячейки гидродинамической сетки, для которой определяются т.н. «эквивалентные» фильтрационноемкостные параметры. Такая процедура в зарубежной литературе называется up’scaling (ремасштабирование) геологической модели в гидродинамическую. В последнее время теория up’scaling’a интенсивно развивается за рубежом. Число отечественных публикаций по этой теме значительно меньше. Завершенная теория up’scaling’a на сегодняшний день еще не построена.

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

Тем не менее, применение даже лучшего из методов upscaling’a не гарантирует от существенных ошибок. Это связано с тем, что погрешность зависит не только от точности расчета эквивалентных параметров, но и от варианта выделения групп мелких ячеек, образующих крупные ячейки.

Процедура выбора вариантов объединения ячеек с контролем погрешности называется upgridding.

Наиболее часто при upscaling’e производят объединение слоев геологической модели. В этом случае сохраняется регулярная структура сетки гидродинамической модели. На практике процедура выбора наилучшего варианта объединения слоев делается вручную, когда из множества вариантов выбирается тот, для которого погрешность, рассчитанная с помощью гидродинамического симулятора, имеет наименьшее значение. Погрешность служит мерой отличия расчетов на геологической и гидродинамической сетках, т.е. до и после upscaling’а. При этом для расчетов выбирается, как правило, не вся геологическая модель месторождения, а специально подобранный «типичный» участок.

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

Целью работы является разработка методик для расчета извлекаемого количества нефти с использованием эквивалентных относительных фазовых проницаемостей при ремасштабировании геолого-гидродинамических моделей пластовых систем и для контролируемого построения гидродинамической сетки.

Научная новизна работы состоит в следующем:

• Предложена оригинальная методика расчета изменения количества извлекаемой нефти из-за ремасштабирования геолого-гидродинамических моделей, основанная на решении задачи двухфазной фильтрации в крупных ячейках, включающая в себя новую процедуру расчета модифицированных ОФП, способ вычисления критических насыщенностей для крупных ячеек гидродинамической сетки, а также формулу для определения коэффициента охвата вытеснением при upscaling’e;

• Разработана новая методика анализа чувствительности геологогидродинамической модели при объединении ее слоев и нахождения наилучшего варианта объединения слоев при upscaling’e, содержащая оригинальный метод расчета погрешности для укрупненных ячеек и процедуру выбора варианта с наименьшей/наибольшей погрешностью.

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

Достоверность результатов работы подтверждается физической непротиворечивостью используемых математических моделей, а также решением тестовых задач, имеющих известные аналитические и численные решения.

Публикации. Основные результаты диссертации опубликованы в работах, список которых приведен в конце автореферата, включая 4 статьи в журналах, рекомендованных ВАК.

Апробация работы. Результаты работы докладывались и обсуждались на следующих конференциях:

1. Международная конференция студентов и аспирантов по фундаментальным наукам «Ломоносов 2005» (Москва, 2005);

2. Международная научно-техническая конференция «Нефть и газ Западной Сибири» (Тюмень, 2005, 2007);

3. 5-ая Всероссийская конференция молодых ученых «Проблемы механики: теория, эксперимент и новые технологии» (Новосибирск, ТатНИПИнефть (Казань, 2006);

5. 50-ая научная конференция МФТИ «Современные проблемы фундаментальных и прикладных наук» (Долгопрудный, 2007) 6. Международная конференция геофизиков и геологов «Тюмень-2007»

(Тюмень, 2007);

7. 52-ая научная конференция МФТИ «Современные проблемы фундаментальных и прикладных наук» (Долгопрудный, 2009);

8. Международная техническая конференция SPE «Повышение нефтеотдачи путем рационального управления заводнением» (Москва, 9. Международная техническая конференция SPE «Взаимодействие геологической и гидродинамической моделей» (Москва, 2010).

Кроме того, результаты работы докладывались и обсуждались на семинарах: в Тюменском государственном университете под руководством профессора А.Б. Шабарова; в Тюменском филиале Института теоретической и прикладной механики им. С.А.Христиановича СО РАН под руководством профессора А.А. Губайдуллина, а также в отраслевых институтах: ОАО «ВНИИНефть», ОАО «СибНИИНП», ООО «КогалымНИПИНефть».

Структура и объем работы. Диссертация состоит из введения, трех глав, заключения и списка литературы. Диссертационная работа изложена на 97 страницах, содержит 29 рисунков и библиографию, насчитывающую наименования.

Благодарности. Автор благодарит своего научного руководителя, д.ф.м.н. С.П. Родионова, а также коллектив Центра моделирования наукоемких технологий нефтедобычи ОАО «СибНИИНП» за сотрудничество.

КРАТКОЕ СОДЕРЖАНИЕ РАБОТЫ

Во введении обоснована актуальность темы диссертации, сформулирована цель работы, отмечена научная новизна и практическая значимость исследований. Диссертация состоит из трех глав.

В первой главе выполнен обзор современных методов upscaling’а и upgridding’а геолого–гидродинамических моделей. Дана классификация методов вычисления эквивалентных параметров для ячеек гидродинамической сетки.

Во второй главе приведено описание численных методов решения уравнений двухфазной фильтрации и выполнено решение тестовых задач.

Изложена предлагаемая процедура вычисления модифицированных относительных фазовых проницаемостей (ОФП). На основе этой процедуры -5исследовано влияние неоднородности распределения фильтрационноемкостных параметров на ОФП и определен коэффициент извлечения нефти (КИН).

В.п.2.1 приведены уравнения двухфазной фильтрации. Изложены численные методы линий тока и IMPES.

В.п.2.2 для проверки правильности работы вычислительной программы для моделирования фильтрационных потоков в нефтяном пласте осуществлено решение тестовых задач, имеющих известные аналитические решения.

В.п.2.3 приведена постановка задачи, метод получения модифицированных ОФП, а также допущения, лежащие в его основе.

Схематическое представление процедуры up’scaling’a для отдельной крупной ячейки приведено на рис.1. Вычисление эквивалентных значений скалярных величин, характеризующих емкостные свойства крупной ячейки, таких как пористость породы и насыщенности фаз Sl не вызывает трудностей и определяется путем осреднения по ее объему V (L.Durlofsky, 2005):

где индексом l внизу отмечены параметры, относящиеся к водной ( l = w ) и нефтяной ( l = o ) фазам. Эквивалентные параметры крупной ячейки отмечены знаком * вверху. Сложнее обстоит дело с определением фильтрационных параметров крупных ячеек, таких как абсолютная проницаемость пористой среды и относительные проницаемости фаз, т.к. они являются не только скалярами, но и векторами или тензорами.

Рис.1. Схематическое представление процедуры up’scaling’a в декартовых координатах x1 и x 2 для крупной ячейки, включающей в себя 88=64 мелкие ячейки; Lk - длины крупной ячейки вдоль координат xk ( k =1, 2).

Рассмотрим, следуя (L.Durlofsky, 2005), процедуру вычисления эквивалентной абсолютной проницаемости в рамках т.н. «локального»

подхода. В каждой точке крупной ячейки используется закон Дарси или в компонентном виде 0 k, вектора скорости жидкости u и градиента давления p ; µ вязкость жидкости. Подставляя выражение (2) в уравнение сохранения массы жидкости получим следующее уравнение для давления p( x1, x2 ) фильтрующейся жидкости, которое решается численно внутри крупной ячейки Тензор эквивалентной абсолютной проницаемости k * для крупной ячейки гидродинамической сетки представляет собой коэффициент пропорциональности между средними по объему скоростью жидкости u и градиентом давления p :

Тензор k * определяется из выражения Здесь скобками обозначены средние величины, являющиеся результатом up’scaling’a, но не используемые в расчетах на крупной сетке.

Знаком * отмечены величины, необходимые для расчетов на этой сетке.

Для решения дифференциального уравнения (3) задаются условия для давления на границах укрупненной ячейки (L.Durlofsky, 2005):

а) – «открытые границы»

Здесь каждому типу граничных условий (а или б) соответствует два варианта их задания (1 и 2).

выбирается тип граничного условия: а) или б). Затем дважды решается уравнение (3) с вариантами j =1 и 2. Система уравнений для определения компонент k * имеет вид:

Здесь индексами 1 и 2 вверху отмечены значения для граничных условий с номерами 1 и 2. Индексами 1 и 2 внизу обозначены соответствующие компоненты вектора p. Чтобы получить симметричный тензор система уравнений (6) дополняется следующим условием:

Из системы (6)-(7) и определяются компоненты симметричного тензора k.

Полный тензор k * может быть определен только при граничных условиях типа а). При условиях типа б) можно вычислить только диагональную часть этого тензора:

Рассмотрим теперь предложенную в диссертации процедуру определения тензоров модифицированных ОФП для водной ( k rw ) и нефтяной ( k ro ) фаз с вязкостями µ w и µo. Внутри крупной ячейки скорости фаз определяются на основе обобщенного закона Дарси (К.С.Басниев, И.Н. Кочина, В.М.Максимов (1993)):

Система уравнений (9) кроме граничных условий (5) дополняется граничными условиями для насыщенности одной из фаз (для воды). Эти условия задаются на границах, через которые происходит втекание жидкости в ячейку. Для вариантов 1 и 2 из (5), имеем:

где S wb (t ) - зависимость водонасыщенности на границе от времени, которая должна выбираться из физических соображений или из приближенных решений. Кроме того, внутри крупной ячейки задаются начальные условия:

Далее, в крупной ячейке рассчитываются средние скорости каждой из фаз Компоненты симметричного тензора модифицированных ОФП l - ой фазы определяются из следующей системы уравнений:

Здесь компоненты тензора k * рассчитываются согласно (6)-(7). Уравнения (12) получены в результате осреднения скоростей фаз, определенных из обобщенного закона Дарси, по объему крупной ячейки.

Полный тензор модифицированных ОФП, как и полный тензор абсолютной проницаемости, может быть определен только для граничных условий типа а). Для граничных условий типа б) значения диагональных компонент усредненных ОФП определяются согласно следующим формулам:

где значения k11 и k 22, определяются из (6) при условии (7).

Алгоритм определения зависимости модифицированных ОФП от водонасыщенности выполняется следующим образом. Внутри крупной ячейки при заданных начальных (11) и граничных (5) условиях на основе системы уравнений (9) решается задача нестационарной двухфазной фильтрации. При этом в каждой точке крупной ячейки и в каждый момент времени рассчитываются насыщенности фаз и давление смеси. Далее, согласно (1) вычисляются средние насыщенности фаз по объему крупной ячейки Sl* = Sl* (t ), средние значения скоростей и градиентов давления, и, согласно (12) или (13), – модифицированные ОФП k rl = k rl (t ). При каждом t значения Sl* и k rl ставятся в соответствие друг другу. В результате получаем искомые зависимости k rl = k rl ( S w ) ( l = w, o ).

Предложенная выше процедура легко обобщается на трехмерный случай.

В этом случае для определения тензоров k * и k rl необходимо решить три задачи двухфазной фильтрации, с тремя компонентами скоростей фаз и градиентов давления и с тремя вариантами граничных условий. Всего, таким образом, определяются девять компонент тензоров абсолютной и относительных фазовых проницаемостей трехмерной крупной ячейки.

В.п.2.4 выполнены расчеты модифицированных ОФП и даны рекомендации для их практического применения на основе процедуры, вытеснением при ремаштабировании.

Рассмотрим некоторые следствия из полученных п.2.3 выражений для модифицированных фазовых проницаемостей (12), (13). Эквивалентная абсолютная и относительные фазовые проницаемости могут быть тензорами, даже если для мелких ячеек они являются скалярами. Это свойство ОФП является следствием неоднородности параметров мелких ячеек. Тензорные модифицированные ОФП не являются функциями только параметров мелких ячеек, а зависят также от параметров жидкостей (например, отношения вязкостей фаз), времени и от конкретной задачи: т.е. начальных и граничных условий. Модифицированные ОФП также зависят от количества мелких ячеек, включенных в крупные ячейки, т.е. от степени укрупнения ячеек.

Поэтому универсального метода up’scaling’a (когда эквивалентные параметры крупных ячеек зависели бы только от параметров мелких ячеек) как для абсолютной проницаемости, так и для относительных проницаемостей не существует. Поскольку эквивалентные параметры крупных ячеек зависят от решения на мелкой сетке, то для получения приемлемого результата необходимо делать предположение об этом решении или использовать приближенные решения.

Обратимся к результатам расчетов эквивалентных ОФП, определенных на основе (13), в некоторых распространенных ситуациях. В качестве примера рассмотрим случай вытеснения нефти водой. Исследование проводилось для двумерной расчетной области в координатах x1 и x2. Для мелких ячеек ОФП задавались согласно модели «разноцветных жидкостей», когда k rl линейно зависят от S rl, а вязкости фаз одинаковы. Численное решение задачи двухфазной фильтрации проводилось методами линий тока и IMPES.

Для давления задавались граничные условия отсутствия потока (тип а), при которых тензор модифицированных ОФП является диагональным. Для водонасыщенности задавались следующие граничные и начальные условия:

Для исходных ОФП S w,min = S wc и S w,max = 1 S or, где S wc и Sor - значения критических насыщенностей, при которых соответствующая фаза теряет подвижность. При граничных условиях (14) средняя водонасыщенность в крупной ячейке S w с течением времени может изменяться во всем подвижном диапазоне от S wc до 1 Sor. В расчетах вязкости фаз полагались равными 1сПз; значения критических насыщенностей S wc и S or для мелких ячеек-коллекторов - 0.2; проницаемость изотропной пористой среды - 100мД;

пористость - 0.2. Для ячеек-неколлекторов эти параметры равны нулю. Во всех ячейках-коллекторах задавался одинаковый вид ОФП.

На рис.2 представлены расчетные зависимости составляющих тензора модифицированных ОФП от средней водонасыщенности крупной ячейки, а проницаемости, имитирующее глинистые включения в песчанике.

Как следует из рис.2, модифицированные ОФП также имеют критические насыщенности S w,min и S w,max. Видно, что значения S w,min и S w,max для ОФП крупных ячеек ближе друг к другу, чем для ОФП мелких ячеек. Это смещение критических насыщенностей обусловлено наличием застойных зон, через которые не проходят линии тока. Величина смещения зависит от основного направления фильтрации.

Рис.2. Зависимости компонент k rl11 и k rl 22 ( l = w, o ) тензоров модифицированных ОФП от средней водонасыщенности крупной ячейки S w (сплошные линии) и ОФП мелких ячеек-коллекторов (штриховые линии). Расчетные точки отмечены знаками, и.

Рис.3. Распределение проницаемости и картины линий тока при расчете компонент тензора krl11 (рис.3а – слева) и k rl 22 (рис.3б – справа). Области коллектора и неколлектора отмечены соответственно черным и белым цветом.

Критические насыщенности для модифицированных ОФП можно получить из следующих соображений. Объем V крупной ячейки равен сумме объемов коллектора и неколлектора. В свою очередь объем коллектора Vc разделяется на подвижный Vcd и неподвижный Vcs. Находящаяся в коллекторе нефть также разделяется на подвижную и неподвижную.

представить в виде:

Здесь первое и второе слагаемые в правой части выражения для *Sl* соответственно подвижная и неподвижная часть порового объема, занятого lой фазой; = Vс / V - коэффициент песчанистости. После несложных преобразований это выражение приводится к следующему виду:

где cd - доля подвижного порового объема коллектора, Sl d и Sl s средние насыщенности l–ой фазы соответственно в подвижной и неподвижной частях коллектора. Величины cd, Sl d и Sl s определяются следующим образом:

соответственно при S w d = S w,min d и S w d = S w,max d. Учитывая, что при постоянных граничных условиях S w s = S w (0) s = S w,min s = S w,max s =const, выражения для S w,min и S w,max принимают вид:

Видно, что S w,min и S w,max зависят от начальных условий. Они зависят также от и граничных условий, поскольку, согласно (17), от них зависит подвижный объем. Однако они не зависят от вязкостей фаз. Для каждого направления значения S wc и S or могут быть различными. Необходимо меньшие значения, чем при «закрытых».

Подставляя (16) и (17) в (18) при начальных и граничных условиях соответственно (14) и S o ( x1, x2,0) = 1 S wc, с учетом того, что S o () = Sor, для критических насыщенностей имеем следующие выражения Отсюда следует, что при отмеченных условиях значения критической насыщенности воды для исходных и модифицированных ОФП одинаковы S w,min = S w,min. Это подтверждается также приведенными на рис.2 численными расчетами.

В плане практического применения описанной выше процедуры можно отметить следующее. В коммерческих программах используется функция водонасыщенности:

т.е. наряду с k rw = k rw ( S w ) можно производить расчеты и с k rw = k rw ( S w ).

Для того, чтобы учесть эффект смещения критических насыщенностей, из тензоров k rl, S wc и S or нужно выделить шаровую часть, а затем применить масштабирование к модифицированным ОФП. Модифицированные нормированные ОФП и водонасыщенность определяются следующим образом:

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

При адаптации гидродинамических моделей к истории разработки месторождения, для которых upscaling уже выполнен, помимо проницаемости можно рекомендовать оценивать и варьировать критические насыщенности S wc и S or. Их нижней оценкой являются величины С помощью критических насыщенностей в диссертации предлагается следующий способ оценки влияния upscaling’a на извлекаемые запасы нефти.

По аналогии с параметрами, характеризующими вытеснение нефти из пласта в целом, введем характеристики степени вытеснения нефти для укрупненных ячеек (или участка пласта). Коэффициент извлечения нефти K по определению равен отношению извлекаемых запасов к начальным запасам в пласте. Структура КИН обычно представляется в виде произведения двух коэффициентов (А.П. Крылов, 1957):

где K выт - коэффициент вытеснения, K охв - коэффициент охвата. Эти коэффициенты определяются соответственно как отношение максимально возможных извлекаемых запасов к начальным запасам и отношение изменения объема нефти в пласте за определенное время к максимально возможному количеству извлекаемой нефти. В свою очередь K охв также представляется в виде произведения двух коэффициентов: коэффициента заводнения K зав и коэффициента охвата вытеснением K охв выт.

Формулу для КИН применительно к крупной ячейке предлагается записать в следующем виде где Здесь S o (0) S o () - максимально возможное количество извлекаемой нефти, получающееся при допущении, что весь объем коллектора крупной ячейки является подвижным. Если S o (0) = 1 S wc, S o () = Sor и S o (0) = S o (0), то S o (0) S o () = 1 S wc S or. Поэтому коэффициент вытеснения для мелкой ячейки выражается через S wc и Sor, а для крупной ячейки - через S wc и S or :

Видно, что S wc, S or в выражении для K выт крупной ячейки играют ту же роль, что и S wc, Sor в аналогичном выражении для мелкой ячейки.

Выражение для K охв выт принимает вид (С.П. Родионов, Л.Н. Соколюк, 2008):

В рамках концепции геолого-гидродинамических моделей определенные выше коэффициенты обладают следующими свойствами. Коэффициент вытеснения является характеристикой внутренних свойств мелких ячеек и не зависит от связей между ними. Коэффициент охвата вытеснением для крупной ячейки напротив можно рассматривать как системную характеристику входящих в нее мелких ячеек. В общем случае, K охв выт, как и модифицированные ОФП, зависит не только от параметров мелких ячеек, но и от начальных и граничных условий, от направления фильтрации. В то же время он не зависит от отношения вязкостей фаз. Это свойство K охв выт позволяет вычислять его по модели «разноцветных жидкостей», которая существенно проще, чем полная двухфазная модель. Поскольку S wc S wc и В третьей главе излагается разработанная в диссертации процедура upgridding’a, в результате которой определяется наилучший вариант объединения слоев гидродинамической модели в процессе upscaling’a, включающая оригинальный способ расчета погрешности. Приводятся некоторые результаты применения upgridding’a.

В.п.3.1 на основе уравнений двухфазной фильтрации рассчитывается погрешность гидродинамической модели при upscaling’e. Используется «локальный» подход, в соответствии с которым погрешность J модели, содержащей Ncell крупных ячеек, представляется в виде суммы погрешностей для отдельных крупных ячеек J = VI (I=1,…,Ncell).

Погрешность I-ой крупной ячейки VI равна сумме (или некоторой норме) x3), т.е. VI = VI, m = ячейки происходит по прямолинейным трубкам тока, проходящим через ряды последовательно расположенных мелких ячеек от поверхностиисточника к поверхности-стоку. Под погрешностью upscaling’а по m-ому направлению (m=13) понимается площадь между зависимостями расхода воды Qw(t) через поверхность-сток крупной ячейки соответственно до Рис.4а-б. Схема потока через крупную ячейку с 123 мелкими ячейками (а) и зависимость для нее расходов воды от времени через границу-сток в направлении m=1 до и после upscaling’а (б). Погрешность VI, m равна площади области, выделенной заливкой на рис.4б.

Погрешность VI, m для крупной ячейки в направлении m=1 может быть выражена через фильтрационно-емкостные параметры входящих в ее состав мелких ячеек. Расход жидкости qi в i–ой трубке тока представляется в виде:

где Tij - проводимость j–ой мелкой ячейки в i–ой трубке в направлении m=1;

Ti - проводимость i–ой трубки; p - разность давлений на нормальных к потоку границах крупной ячейки; k1,ij, s1,ij, x1ij - для ij-ой ячейки соответственно проницаемость, площадь поперечного сечения и размер вдоль декартовой координаты x1. Объемы подвижного порового пространства ij-ой мелкой ячейки ( V p,ij ) и I–ой крупной ячейки ( V p, I ) связаны с расходом жидкости qi и временами движения через них фронта вытеснения в следующем виде: V p, ij = qi ij, V p, I = QI I, где Здесь V p,ij и V p, I выражаются через геометрический объем Vij, пористость ij и критические насыщенности S wc,ij, S or,ij следующим образом:

поровый объем, расход жидкости, времена прохода через ячейки фронта вытеснения проводимости и разности давлений, выражение для погрешности крупной ячейки преобразуется к виду (С.П. Родионов, Л.Н.Соколюк, 2009):

Остальные компоненты погрешности VI вычисляются аналогично VI,1, определяемой из (19).

В отличие от работ (M.King et al, 2006) и (S. Hosseini et al., 2008), вычисление погрешности в диссертации основано непосредственно на решениях уравнений двухфазной фильтрации.

Погрешность гидродинамической модели может быть минимизирована на множестве вариантов объединения мелких ячеек в крупные (или для выбора из определенного множества гидродинамических сеток сетки с минимальной погрешностью).

В.п.3.2 приводится описание разработанного метода для поиска наилучшего (с точки зрения наименьшей погрешности) варианта объединения слоев геологической модели при upgridding’е, содержание которого заключается в следующем.

Пусть геологическая модель (т.е. модель до upscaling’a) состоит из N последовательно прилегающих друг к другу горизонтальных слоев ячеек с номерами l=1,2,..,N0. Слои в декартовой системе координат расположены по вертикальной оси x3, причем каждый слой содержит по оси x3 одну ячейку.

Количество ячеек по каждой из горизонтальных осей x1 и x2 одинаково для всех слоев. Множество вариантов объединения N0 слоев геологической модели в N (1NN0) слоев гидродинамической модели характеризуется матрицами W N, состоящими из N строк и 3-х столбцов:

Здесь и далее индексом N вверху отмечены величины, относящиеся к гидродинамической модели, содержащей N слоев. Строка с номером l матрицы W N, Wl N = l wlN2 wlN3, означает, что l-ый слой (1lN) гидродинамической модели включает в себя с wlN2 -го по wlN3 -ый слои геологической модели. Элементы матриц W обладают следующими свойствами:

Первый элемент второго столбца и последний элемент третьего столбца при любом N принимают значения соответственно 1 и N0. Поскольку матриц W определяется натуральными числами wlN3 (l=1,…,N), из которых N-1 являются независимыми. Матрица W соответствует преобразованию геологической модели самой в себя, а матрица W 1 - преобразованию геологической модели в гидродинамическую модель с одним слоем, т.е. в двумерную модель. Эти две матрицы известны априори.

Погрешность модели в целом представляется в виде целевого функционала, равного сумме погрешностей (функционалов) по отдельным крупным ячейкам. Поэтому функционал модели, состоящей из N слоев ( J N ), равен сумме функционалов по отдельным слоям, а функционал для отдельного слоя ( J lN ), в свою очередь, равен сумме функционалов для входящих в него крупных ячеек:

Здесь VI - погрешность для крупной ячейки с номером I, образовавшейся в результате объединения с wlN2 -го по wlN3 -ый, слоев геологической модели в l–ый слой гидродинамической модели. Из (21) следует, что функционал N-слойной гидродинамической модели определяется матрицей W N ( J N = J N (W N ) ), а функционал l–го слоя этой модели - l–ой строкой Wl N матрицы W N ( J lN = J lN (Wl N ) ).

При объединении слоев ячейки укрупняются только по вертикали, в направлении 3 (рис.5). Количество мелких ячеек в крупной ячейке в направлениях 1 и 2 в этом случае обычно равно n1=n2=1, а в направлении 3 n3>1.

Рис.5. Схематическое представление одного из вариантов объединения ячеек геологической модели, состоящей из 5-ти слоев, в ячейки гидродинамической модели из 3-х слоев. Количество мелких ячеек в слое из крупных ячеек по направлениям 13:

n1=n2=1, n3=3.

Формула для VI в трехмерном случае применительно к рассмотренному выше сценарию объединения слоев с учетом (19) принимает вид:

где Здесь для i3-ой мелкой ячейки: Tm,i 3, k m,i 3, sm,i 3, xm,i 3 - соответственно проводимость, проницаемость, площадь поперечного сечения и размер вдоль координаты xm (m=1, 2); V p,i 3 - подвижный объем; Vi 3 - геометрический объем;

i 3 - пористость; S wc,i 3, S or,i 3 - критические насыщенности фаз.

С помощью матриц W N можно оценивать качество up’scaling’a. Для этого необходимо подставить в (21) предполагаемые варианты объединения слоев и сравнить значения соответствующих функционалов. Однако наибольший интерес представляет определение такой матрицы W N, при которой функционал (21) принимает минимальное значение на определенном множестве объединений слоев. В диссертации предлагается следующая процедура решения этой задачи.

Объединение N0 слоев в N слоев можно представить в виде последовательного выполнения i=1,…,N0-N шагов. На i–ом шаге объединяются только какие-либо два соседних слоя, которые выбираются из N(i-1)=N(i)+1 слоев, имевшихся на предыдущем i-1–ом шаге. Поскольку на данном шаге число слоев на один меньше, чем на предыдущем, то на i–ом шаге имеется N(i)=N0-i слоев (i=1,…,N0-N) и N0-i-1 возможных вариантов объединений двух соседних слоев.

Если вариант объединения двух слоев в один выбран, тогда остальные слои не объединяются и остаются такими же, как и на предыдущем шаге.

Следовательно, соответствие номеров слоев на i-1–ом и i–ом шагах, когда слои l(i-1) и l(i-1)+1 объединяются в один l*(i)-ый слой, можно представить в виде где номера l(i) и l(i-1) могут принимать значения из интервалов соответственно l(i)=1,…,N(i) и l(i-1)=1,…, N(i-1). Поэтому когда на i -ом шаге объединяются два слоя с номерами l(i-1) и l(i-1)+1 в один слой с номером l*(i), строки матрицы W N (i ), неизвестной на данном шаге, выражаются через элементы матрицы W N (i1), известной на предыдущем шаге, следующим образом:

Видно, что при переходе от i-1-го к i-му шагу строки матриц и значения функционалов для всех слоев за исключением l*(i) –го не изменяются. С учетом этого обстоятельства и выражения (21) можно показать, что приращение функционала для модели в целом из-за образования нового l*(i)-го слоя, получающегося в результате объединения l(i-1)-го и l(i-1)+1-го слоев, равно:

Если в качестве управляющей переменной, определяющей тот или иной вариант объединения слоев на i–ом шаге, принять порядковый номер слоя l=l*(i), то вектор L*=(l*(1),…,l*(i),…,l*(N0-N)), составленный из значений l*(i) (i=1,…,N0-N) будет определять вариант объединения N0 слоев в N слоев.

При известном векторе L* матрицы W N и функционалы J N (N=N0с помощью (20) и (21) могут быть легко определены пошагово.

Отметим, что неизвестные (или свободные) компоненты вектора L*, могут N)=l*(N0-1)=1.

Если на каждом шаге i=1,…,N0-N номера объединяемых слоев L* выбирать, в частности, таким образом, чтобы элемент вектора приращений функционала J lN (i ) при l=l*(i) принимал наименьшее значение из множества приращений J lN (i ) (l=1,…, N(i)), т.е:

то функционал при объединении N0 слоев в N слоев будет принимать минимальное значение. Множество вариантов выбора пар слоев для объединения на каждом шаге определяется следующей матрицей:

Приращение функционала J lN (i ) = J lN (i ) (Wl N (i 1) ) для каждого из этих вариантов (l=1,…,N(i)) определяется выражением (23).

Для определения наилучшего варианта объединения N0 слоев в N слоев предлагается использовать следующий «жадный» алгоритм, состоящий из N0-N шагов (N0>2, N>1).

и начальные значения погрешностей J l = 0 (l=1,…,N0). В соответствии с (24) рассчитывается вектор W N 1 для (N0-1) возможных вариантов объединения рассчитываются приращения функционала J lN (i ).

Шаг i=i+1. Число слоев уменьшается на единицу (NN-1). Из множества l=1,…,N(i) определяется номер слоя l=l*(i) при котором J lN (i ) как функция l принимает наименьшее значение. Слои l(i-1) и l(i-1)+ объединяются и новому слою присваивается номер l*(i). Рассчитывается новый вектор wlN3(i ). Для вновь образовавшегося слоя рассчитывается компонента вектора J lN((ii)). Если N>1, то для двух или одного вариантов объединения слоев, граничащих с этим новым слоем, согласно (24) и (23) пересчитываются значения Wl N (i ) и J lN (i ). Осуществляется переход к следующему шагу.

В.п.3.3 приводится демонстрация применения разработанного в п. 3.2.

метода по объединению слоев. На рис.6 приведено схематическое представление работы алгоритма для геологической модели, содержащей четыре слоя. На первом шаге объединяются слои 2 и 3, на втором шаге – слои 1 и 2 и на третьем шаге – 1 и 2. Заливкой отмечены слои, образующиеся в результате объединения слоев с предыдущего шага.

Рис.6. Схематическое представление работа алгоритма для модели из четырех слоев Для случая рис.6 начальная матрица W4 и матрицы W соответственно на первом, втором и третьем шагах выглядят следующим образом:

Вектор L* имеет вид: L*=(l*(1), l*(2)), при этом l*(1)=2, l*(2)=1. Число слоев на каждом шаге принимает значения N(0)=4, N(1)=3, N(2)=2, N(3)=1.

Если число слоев N не является фиксированным, то предлагаемая процедура может играть вспомогательную роль для определения необходимой детальности гидродинамической модели. Поскольку на каждом шаге число слоев уменьшается на 1, то в процессе выполнения шагов имеем последовательность гидродинамических моделей от исходной (геологической), с числом слоев N=N0, до двумерной при N=1. Если сравнить минимальные значения функционалов для каждой из этих моделей то по виду функции J ( N ) = J min можно выбрать необходимую детальность гидродинамической модели. Оптимальное число слоев N выбирается исходя из компромисса между временными затратами на проведение расчетов и их точностью. В качестве безразмерного критерия погрешности можно использовать долю порового объема всего пласта J min / V p или относительную погрешность J ( N ) = J min / J 1 (0 J (N ) 1 при 1NN0).

погрешностей от числа слоев N для геологической модели, состоящей из N0=14 слоев (б), имеющей профиль абсолютной проницаемости k (мд) 0.00864 (а). Точки, отмеченные на рис.7б цифрами 13 соответствуют объединению пар слоев на рис.7а.

Это демонстрируется на следующем примере upscaling’a геологической модели из 14-ти слоев. В ячейках каждого слоя задается постоянное значение проницаемости, но это значение для разных слоев – разное (см. рис.7б).

Размеры ячеек и пористость во всем пласте одинаковы.

На рис.7б представлена зависимость J min / J 1 от N. Видно, что слои, находящиеся в области экстремумов изображенной на рис.7а функции распределения относительной проницаемости, можно объединить без значительной потери точности. Для оценки чувствительности решения и для выявления слоев, наиболее неподходящих для объединения, можно использовать зависимость J max / J 1, также изображенной на рис.7б.

ОСНОВНЫЕ РЕЗУЛЬТАТЫ И ВЫВОДЫ РАБОТЫ

1. На основе развитой процедуры определения модифицированных относительных фазовых проницаемостей (ОФП) установлено, что в неоднородном пласте при объединении мелких ячеек геологической сетки в крупные ячейки гидродинамической сетки (up’scaling), модифицированные ОФП являются тензорными функциями, даже если ОФП мелких ячеек – скалярные функции. Модифицированные ОФП параметров входящих в их состав мелких ячеек, но и от начального распределения насыщенностей фаз и направления фильтрационного потока в крупных ячейках, а также от условий на их границах. В частности, значения критических насыщенностей при «открытых»

граничных условиях меньше, чем при «закрытых».

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

3. Разработана оригинальная методика контролируемого построения гидродинамической сетки, которая позволяет: а) оценивать степень пригодности вариантов для объединения слоев геологической модели в слои гидродинамической модели с учетом погрешности огрубления; б) решать задачу минимизации погрешности гидродинамической модели на множестве вариантов объединения мелких ячеек в крупные при заданном числе слоев; в) получать зависимость минимальной/максимальной погрешности от числа слоев, из которой можно определять необходимую детальность гидродинамической модели, а также оценивать диапазон варьирования числа слоев для процедуры up’scaling. На конкретном примере показано, что с помощью оптимизации процедуры up’scaling можно без потери точности сократить до 50% слоев геологической модели.

1. Родионов С.П., Орехова Л.Н. Расчет относительных фазовых проницаемостей при преобразовании геологической модели в гидродинамическую // Тезисы 50-й Юбилейной научной конференции МФТИ, Долгопрудный – 2007, С.97-98.

2. Родионов С.П., Орехова Л.Н. Расчет относительных фазовых проницаемостей при преобразовании геологической модели в гидродинамическую // Тезисы Всероссийской научно-технической конференции «Нефть и газ Западной Сибири», Тюмень. – 2007, С.45-47.

3. Родионов С.П., Орехова Л.Н. Определение модифицированных относительных фазовых проницаемостей при преобразовании и газ. – 2008. – №6. – С.12-17.

4. Родионов С.П., Орехова Л.Н. Определение модифицированных относительных фазовых проницаемостей при преобразовании геологической модели в гидродинамическую. Ч.2 // Известия вузов. Нефть и газ. – 2009. – №1. – С.4-9.

5. Родионов С.П., Орехова Л.Н. Методика расчета фильтрационных характеристик при апскейлинге геологических моделей. // Вестник ТюмГУ. – 2008. – №6. – С.64-69.

6. Родионов С.П., Соколюк Л.Н. Методика контролируемого up’scaling’a геологических моделей, основанная на упрощенном решении уравнений двухфазной фильтрации // Труды 52-ой научной конференции МФТИ, Долгопрудный. – 2009. – Т.2, Ч.3. – С.60-61.

7. Родионов С.П., Соколюк Л.Н. Расчет и использование модифицированных относительных фазовых проницаемостей при преобразовании геологической модели в гидродинамическую // Труды МФТИ. – 2010. – Том 2, №2. – С.130-136.

1. Louis J. Durlofsky. Upscaling and Gridding of Fine Scale Geologocal Models for Flow Simulation. Department of Petroleum Engineering, Stanford University, Stanford, CA 94305-22220 USA, 2005.

2. Басниев К. С., Кочина И. Н., Максимов В. М. Подземная гидромеханика:

Учеб. – М.: Недра, 1993.

3. King M.J., Burn K.S., Wang P., Venkataramanan, M., Alvarado, F., Ma, X., and Datta-Gupta A. SPE Reserv. Eval. Eng., 9, 2006, 317-334.

4. Hosseini S.A. and Kelkar M. Paper SPE 116113.