Математическое моделирование и параллельные вычисления процессов тепломассопереноса при эксплуатации единого комплекса «Нефтяной пласт - система скважин - электроцентробежные насосы»

УДК 519.6

В. М. Конюхов, И. В. Конюхов, А. Н. Чекалин

МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ И ПАРАЛЛЕЛЬНЫЕ ВЫЧИСЛЕНИЯ

ПРОЦЕССОВ ТЕПЛОМАССОПЕРЕНОСА ПРИ ЭКСПЛУАТАЦИИ ЕДИНОГО КОМПЛЕКСА

«НЕФТЯНОЙ ПЛАСТ - СИСТЕМА СКВАЖИН - ЭЛЕКТРОЦЕНТРОБЕЖНЫЕ НАСОСЫ»

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

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

Keywords: numerical simulation, filtration, oil reservoir, well system, heat- and mass transfer, multiphase flow, centrifugal

electric pump.

An united mathematical model of interconnected heat and mass transfer in the non-uniform porous oil reservoir covered by arbitrary located injecting and producing wells equipped with electrical submersible pumping systems is developed. The parallel calculation technologies are used in numerical computations and visualization of its results. An availability of the united model and significant decreasing of the problem time at parallel computing in comparison with sequential ones are demonstrated by examples of the numerical experiments.

Введение

Нефтяное месторождение и система вскрывающих его нагнетательных и добывающих скважин представляют собой единый объект, работа всех элементов которого тесно взаимосвязана и сопровождается сложными процессами тепломассопереноса. Нефтяные пласты, как правило, имеют сложное геометрическое строение и неоднородную по абсолютной проницаемости и пористости структуру. Процесс фильтрации в пласте представляет собой течение многофазных многокомпонентных смесей, в котором происходит нетривиальное взаимодействие фаз, а нефти могут проявлять неньютоновские свойства. В свою очередь, подъем пластовой продукции на поверхность земли с помощью добывающих скважин также сопровождается нелинейными термо- и гидродинамическими эффектами, фазовыми переходами, изменением структурных форм течения смеси и т. д. как в трубах скважин, так и в каналах подземного насосного оборудования, например, установок погружных электроцентробежных насосов (УЭЦН), применяемых для перекачки водонефтегазовых смесей. В состав таких установок входят ЭЦН, погружные электродвигатели (ПЭД), газосепараторы и диспергаторы, работа которых также зависит от состава перекачиваемой продукции, расходов фаз и других температурных и гидродинамических факторов.

Для управления и диагностики работы подземного оборудования (например, ПЭД) широко используются наземные станции управления (СУ). Они имеют не только прямую связь с УЭЦН, инициируя включение или отключение двигателя установки, но и осуществляют прием и обработку энергетических параметров (силы тока ПЭД, сопротивления жилы питающего установку кабеля и пр.) и показаний телеметрических датчиков УЭЦН о температуре и давлении на приеме ЭЦН, степени виб-

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

Математические модели, применяемые для описания процессов тепломассопереноса в отдельных элементах комплекса «пласт - система скважин - УЭЦН» (см., обзоры в работах [1-6]), строятся с учетом конкретных специфических особенностей многофазных потоков в областях их движения -соответственно, в пористой среде, трубах скважины, проточных каналах ступеней насосного узла погружной установки. Эти модели, как правило, сами по себе представляют собой сложные системы нелинейных дифференциальных уравнений, которые существенно отличаются как по типу, так и по методам их решения. Исследование же режимов эксплуатации всего комплекса как единого целого и их оптимизация требуют совместного решения этих систем уравнений, что приводит к значительным трудностям. Поэтому большинство публикаций (например, [7,8,9-13]) по данной тематике посвящено вопросам изучения характеристик отдельных элементов комплекса без учета их взаимосвязи. Так, при расчете фильтрационного процесса в нефтяном пласте учет работы добывающих и нагнетательных скважин, вскрывающих коллектор, моделируется граничными условиями, в которых на забое каждой скважины задается либо давление, либо ее суммарный дебит. В первом случае в результате решения фильтрационной задачи определяется суммарный дебит, а во втором - забойное давление. При этом совершенно не учитывается тот факт, что условия на забое могут оказаться недостаточными или вообще неприемлемыми для обеспечения нормальной работы скважин.

Такие ситуации могут возникать, например, в добывающих скважинах, в которых для подъема продукции используются погружные электроцентробежные установки. Это связано с тем, что допустимый диапазон подач, в котором работа ЭЦН является стабильной, достаточно узок: их значения должны составлять примерно от 0,7 до 1,3 долей от оптимальной производительности насоса при перекачке многофазной смеси. В свою очередь, рабочий диапазон подач и напорная характеристика ЭЦН существенно зависят от структуры и состава смеси на приеме насосного узла. Состав смеси определяется как обводненностью продукции, поступающей в скважину из пласта, так и процессами межфазного массообмена при ее движении в лифтовых трубах. Массообменные процессы возникают по мере снижения давления в скважине ниже давления насыщения нефти газом, приводя к образованию новой свободной газовой фазы в потоке. Если величина забойного давления такова, что разгазирование происходит ниже координаты подвески насоса, то на его прием поступает смесь, наличие газа в которой при объемной газонасыщенности свыше 25% на входе в первую ступень оказывает весьма вредное влияние на характеристики работы ЭЦН (вплоть до срыва подачи при газосодержании ~35%).

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

Сложность численного моделирования процессов в едином комплексе «пласт - система скважин - УЭЦН» обусловлена не только специфическими особенностями процессов в каждом из ее элементов, но и затратами машинного времени для проведения расчетов. Однако в последние годы происходит бурное развитие многопроцессорных компьютеров и суперЭВМ, технологий параллельных вычислений и алгоритмических языков программирования. Их применение позволяет повысить эффективность реализации подобных задач на ЭВМ и связанных с ними задач анализа и прогнозирования рабочих характеристик всех элементов системы разработки нефтяного месторождения.

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

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

Особенность построения общей математической модели заключается в необходимости описания процессов, происходящих при течении многофазных смесей в различных областях комплекса «пласт - система скважин - УЭЦН», и последующего сопряжения решений соответствующих уравнений на границах этих областей.

Математическая постановка

Математическая модель фильтрации в пласте в рамках крупномасштабного приближения. Рассмотрим процесс двухфазной двухкомпонентной (вода, нефть) фильтрации в неоднородном нефтяном пласте, достаточно тонком по сравнению с его размерами в плане, вскрытом системой произвольно расположенных скважин. Абсолютная проницаемость К, пористость m, толщина Н, начальная водо-насыщенность 50 могут изменяться по простиранию коллектора, его кровля и подошва непроницаемы. Так как мощность пласта значительно меньше его протяженности, то осредненная по толщине коллектора система нелинейных дифференциальных уравнений [1,2] фильтрационного процесса в области 0пп=(0<х<1х, 0<х<Ьч} при 0<<ТР без учета капиллярных эффектов, силы тяжести, сжимаемости фаз и пористой среды в изотермическом приближении может быть записана в виде

- = 0,

сV + д_У1_

дх ду

дР

V = - ккф н—,

дх

V = -НККф

дР

ду

дХ (^)+-

дх

д5

(V) + тН— = 0 =

у дг

КФ = КФ1

КФ1 =

ду

/М + кфз/Мз , / = Кфз/(/"гКф) ,

(5* - 5): /(5* - ) ,0 < 5 < 5*,

0 , 5 * < 5 < 1,

кфз

(5 - 5* )3/ 5:

0 < 5 < 5*, 5* < 5 < 1,

где г, х и у - временная и пространственные координаты; Ьх и Ьу - протяженность области Бпл вдоль осей Ох и Оу, Тр - время разработки пласта; Н=Н(х,у), т=т(х,у), К=К(х,у);

\(х,у,г)=\1(х,у,г)+\2(х,у,г), 5(х,у,г) и Р(х,у,г) - суммарная скорость фильтрации, водонасыщенность и давление; У|, Кф\, м - скорость, относительная проницаемость и вязкость 1-ой фазы (здесь и в дальнейшем нижние индексы «_/», «2» и «з» обозначают параметры нефти, газа и воды соответственно); /доля воды в суммарном фильтрационном потоке; 5* и 5 - значения связанной и предельной водонасы-щенности; :1 и :3 - показатели степенных зависи-

мостей фазовых проницаемостей. В общем случае

«1 ^ аз.

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

Трехфазный поток в трубах добывающей скважины. Квазистационарные процессы тепломас-сопереноса при движении газоводонефтяной смеси в скважине могут быть описаны системой уравнений [3-5], полученной в результате осреднения и конкретизации общих моделей [6] с учетом специфических особенностей скважинных потоков. Эта система включает в себя законы сохранения импульса и энергии трехфазного потока в обсадной колонне 0ок—{0<х<1э} и насосно-компрессорных трубах (НКТ) Б нш—{2э<2<Ис}, а также дополнительные соотношения, связывающие между собой различные характеристики потока:

дР 2т,

дх

д 3

— У Ом= / дх 1=1 ' 1

(

\

----— + pg  008 ф

У ОС

дТ т^ О р ¡"Г" = Т У ар 1

дх ^Р

дР_ ~дг'

+2лГт (тСтъ - дст) + /т (qv - 3Ь) >

3

p = Уpm, ™ = сi+vi, ] = 2,3,

1=1

РплРдег

Р1 =

Сн =

Рз =

Р дег

Р2,(Л Р

Р1п

Р2 =

^г ят

Рз =Рзп

(Р,Т),

Рг =•

Оз (Р (1-&2 ) + Р2&) , Рз (О - О3 ) + РО3 ' Щ _О2Р1Рз_

Рз [Р02 + Р2 (О - 02 )] + ОзР2 (Р - Рз ) Здесь х - координата пространственной вертикальной оси Ох с началом на кровле пласта; Ис - глубина скважины; р, О', щ и Р - средние по сечению / труб скважины плотность, средневзвешенная скорость, массовый расход, истинное и расходное объемные содержания / -ой фазы; Р, Т, р, V, О—О1+О2+Оз - давление (одинаковое для всех фаз), температура, плотность, объемная скорость и массовый расход смеси; ар, ат - коэффициенты теплового расширения и объемной упругости /-ой фазы; Я - универсальная газовая постоянная, 2г - коэффициент сверхсжимаемости реального газа; £(Р,Т)—1-аР' (Т-Тпп) + ап (Р-Рпп);Рпп и Тпп - начальные пластовые давление и температура; р1пп и рзпп - значения плотности нефти и воды в пластовых (Р—Рпп, Т—Тпп); р0 и р-|дег - плотность газа и дегазированной нефти в нормальных условиях (при атмосферном давлении Рат и температуре Т20 —20 °С; дф) -

отнесенная к единице объема смеси линейная мощность распределенных вдоль НКТ внешних источников тепла, обусловленных потерями электроэнергии в кабеле, питающем УЭЦН; гст и дст - средние по периметру канала касательное напряжение и плотность теплового потока на стенке скважины; 3 -массовая скорость образования газовой фазы за счет нефтяной, отнесенная к единице объема смеси; Ь -скрытая теплота растворения газа в нефти; Сн - массовая концентрация газа, растворенного в нефти при давлении Р выше давления Рн насыщения нефти газом; Я(Р/Рн,Т) - обобщенный коэффициент растворимости газа в нефти; У0(Т) - нормальный объем газа, выделяющийся из нефтяной фазы при Р=0; С] -аналоги коэффициента Зубера, учитывающие неравномерность распределения по сечению /т скорости смеси и объемных истинных концентраций ]-ой дисперсной фазы; ^ - ее скорость дрейфа; гт - радиус канала (значения т=«0к» и т=«нкт» соответствуют участку Б0к обсадной колонны скважины ниже точки 2э подвески погружной установки и насосно-компрессорным трубам Бнкт); g - проекция вектора ускорения силы тяжести на ось Ох; ф(х) - угол наклона профиля скважины относительно оси Ох.

Математическая модель (6)-(8) учитывает такие основные эффекты в трубах скважины, как неизотермичность, фазовые переходы, трение, влияние силы тяжести, разгазирование нефтяной фазы в скважине, сжимаемость фаз, изменение структурных форм течения смеси, возможную инверсию жидких фаз и проскальзывание дисперсных составляющих потока, его теплообмен с окружающей средой.

Трехфазный поток в каналах ЭЦН. Для расчета характеристик процесса тепломассопереноса при движении водонефтегазовой смеси в насосном узле Бэцн—{0<%<Ьэ} и параметров работы отдельных его ступеней воспользуемся математической моделью [4,5,9], обобщающей дискретную схему П. Д. Ляпкова и В.И. Игревского [10]:

, дР 1 з дт

1с — = gp^И, — У О'Срг—=

/э и

V 1=1 Р V Здесь % - координатная ось с началом на входе в первую ступень, проходящая вдоль оси насосного узла; Ьэ=1С  Мс - общая длина пакета ступеней, 1С -длина отдельной ступени, Мс - общее число ступеней в насосном узле; /э=п гэ2, гэ и иС=/э1С - эффективные радиус каналов и объем проточной области ступени; ЛИ, г]=gpЛИQ/ЛN, Q—G/p и ЛЫ - напор, коэффициент полезного действия, объемный расход и потребляемая мощность, характеризующие работу отдельной ступени при перекачке смеси; KQ=Q /рв, Кн—ЛИ/ЛИв и К—ц/Цв. - коэффициенты относительных подачи, напора и кпд. работы ступени, отражающие снижение ее параметров Q, ЛИ, V при перекачке вязкой среды по сравнению с соответствующими показателями Qв, ЛИв, гв работы ступени на воде; д—Q/Qв - относительный расход смеси.

- Г т У Р

^рР^лИ -3Ь,

Для определения параметров ЛИ, г и ЛИ при работе ступени ЭЦН на неоднородной смеси воспользуемся расчетными формулами, полученными в работе [11] на основе модификации методики П. Д. Ляпкова [12]. Так, в диапазоне относительных подач 0.5 < дв = Qв/< 1.3 для ламинарного и для

турбулентного режимов течения напор ЛИ ступени при перекачке смеси расходом < и эффективной вязкостью м с погрешностью 1 ^ 2% вычисляется по следующим формулам:

п, = 193^72* (ли;у75,

/(и8) = (4.3 + 0.816п0274) п50 575 , г = м1р, Ке = /(п, Ка = Кн = 2(К6- 200*У ,

2в = 2Ка , ЛИ = КнЛИв

Здесь 2*, лИ*, г* - паспортные значения подачи, напора и кпд при работе ступени ЭЦН в оптимальном режиме на воде при #в=1; <в, ЛНВ, /в - аналогичные характеристики при дв ^ 1, определяемые по паспортной кривой «напор - подача» ступени, п5 -коэффициент ее быстроходности, учитывающий конструктивные особенности проточных каналов ступени; ю - угловая скорость вращения вала насоса (ротора электродвигателя), V - коэффициент эффективной кинематической вязкости смеси. Величина м зависит от концентрации дисперсных фаз в потоке и определяется по формулам [4,5] с учетом эффекта инверсии фаз. Величину мощности ЛИ, потребляемой ступенью ЭЦН, и ее кпд г при перекачке эмульсии в диапазоне 0.5<дв<1.3 определим соотношениями [12]:

лМ = лыв

Кеэ = Ка Ке 1в , Кг = т1п (1,т1п (Кг,л, Кг,Т )),

Кг,т = 0.2741В Кеэ - 0.06 - 0.Ща,

Кг,л = 0.485^ Яеэ - 0.63 - 0.26^,

где Кцл, КгТ -коэффициенты относительного кпд ступени при ламинарном и турбулентном режимах течения смеси, Яеэ - число Рейнольдса, введенное П. Д. Ляпковым. Зависимости ЛИ(2в) и ЛИв(<в) задаются как табличные функции подачи <в.

Для пересчета характеристик ступени с воды на неоднородную смесь при 0<дв<0.5 и 0.5<дв<двтах построены зависимости [9], аналогичные (11), (12), где двтах - относительная максимальная подача ступени по воде, при которой напор ЛИв равен нулю.

Возможность изменения режима работы УЭЦН с помощью станции управления (например, «ЭЛЕКТОН») учтем при помощи зависимости параметров работы насоса и мощности МД, потребляемой погружным электродвигателем, от частоты [14]:

2в = 2вн /н , ЛИ в =ЛИвН

ЛМв = ЛМвн (Ю/Юн )3, Мд = Мдн /н ,

где <вн, ЛИвн, гвн и Мвн - характеристики ступени и двигателя установки при номинальной угловой ско-

рости вращения вала ПЭД (юн=50 Гц), а <в, ЛИв, гв и Мд - при частоте со, отличной от юн.

Замечание 1. Модели (6)-(8) и (9), (10) как частные случаи включают в себя уравнения однофазного (нефтяного или водного) и двухфазного (газонефтяного или водонефтяного) потоков в трубах скважины и каналах ЭЦН. Например, уравнения течения газонефтяной смеси получаются из уравнений (6)-(8) и (9), (10) при 03=0.

Замечание 2. В данной статье приведены лишь некоторые основные соотношения для характеристик потоков, рабочих параметров ступеней и погружного двигателя. Полное описание всех полуэмпирических зависимостей, аппроксимационных формул и расчетных соотношений, необходимых для замыкания моделей (6)-(8) и (9)-(13) можно найти в работах [3-5, 9, 11, 13].

Замечание 3. Уравнения (6)-(8) и (9)-(13) записаны для отдельной добывающей скважины _Оскв, содержащей конкретную насосную установку _Оэцн. В рассматриваемом случае нефтяной пласт вскрыт системой добывающих и нагнетательных скважин, количество которых равно соответственно Мн и Мд. Очевидно, что полная математическая модель содержит Мд таких систем уравнений для расчета характеристик работы каждой к-ой добывающий скважины (к=1,..Мд), отличающихся друг от друга своими геометрическими параметрами (инк-линометрическим профилем, глубиной Ис,к, радиусами гок,к и гнктк труб обсадной колонны и НКТ), типами установок электроцентробежных насосов и погружных электродвигателей, координатами спуска УЭЦН в скважину, составом поступающей из пласта продукции, распределениями термо- и гидродинамических параметров, величинами Рз,д,к и Ру,д,к забойных и устьевых давлений и т.д. Для упрощения записи в уравнениях (6)-(8), (9)-(13) индекс «к», обозначающий порядковый номер конкретной скважины, опущен.

Замечание 4. Для расчета давления Рзн на забое каждой ]-ой нагнетательной скважины 0=1,..Мн), вообще говоря, следует решать Мн систем уравнений, аналогичных уравнениям (5)-(6) и описывающих течение водяного потока при заданной величине РН,уст давления закачки на устье таких скважин. Однако, как показано в работах [3,5], в данной ситуации можно пренебречь сжимаемостью и тепловым расширением воды, а также силами инерции и трения. Тогда в рамках изотермического приближения величина Р3Н определятся из простого уравнения гидростатики, так что Рз,н=Ру,н+ярИс, где Ру,н - заданное давление на устье >ой нагнетательной скважины (здесь также опущено индексное обозначение параметров нагнетательных скважин).

Сопряжение математических моделей. Решение фильтрационной задачи (1)-(5) отыскивается при 0<г<Тр в прямоугольной области Бил, в которой в начальный момент времени г=0 задано распределение водонасыщенности:

5

Поскольку уравнение (3) относительно насыщенности является уравнением переноса, то эта функция

должна быть задана на боковых поверхностях Гн, нагнетательных скважин и на той части границы Г, через которую поток жидкости направлен внутрь области решения. Обозначим эту часть границы через Г(+) и запишем соответствующие граничные условия:

5

= 5 (+>, (5* < 5(+) < 5 *),

5 ( X, У )ЕГн,= 5

} = 1, М н.

Уравнения (1), (2) служат для определения давления в пласте. Нетрудно видеть, что подстановка (1) в (2) приводит к уравнению эллиптического типа относительно функции Р, для которого на внешней границе Г области Бпп задается давление Рг:

Нагнетательные и добывающие скважины с границами Гн,|, Гдк радиусом гок,|, гокк работают в режиме заданных забойных давлений. Соответствующие граничные условия имеют вид:

Р

= Р

(х.у)еГн, з,н,]

] = 1, М н:

Р\ = Р

(х,У)еГлк з,д,к '

к = 1, М„

Если известно начальное распределение насыщенности 5о(х,уЗ) в пласте и заданы величины Рг, Рз,н,1, Рз,д,к, 5(+), 5», 5, то фильтрационная задача (1)-(5), (14)-(17) может быть решена независимо от системы уравнений в скважине. В результате будут найдены поля давления Р(х,уЗ), водонасыщенности 5(х,уЗ), суммарный поток У(х,уЗ) и потоки У 1(х,у,1), У2(х,у^) отдельных фаз, а также такие важнейшие характеристики работы скважин, как величины суммарного объемного расхода Р двухфазной смеси и объемного расхода р1 воды через поверхность Гд вдоль направления внешней нормали п, а также обводненность в дебита добывающей скважины:

д . ^ = 0 /ккАу, в = Р. (18) дп ГД дп р1

При моделировании процессов тепломассо-переноса в едином комплексе «пласт - система скважин» для сопряжения решений систем уравнений (1)-(5) и (6)-(8) должны быть заданы условия, связывающие характеристики фильтрационного потока на боковой поверхности Гдк забоя каждой кой добывающей скважины и потока в обсадной колонне Бок на уровне кровли пласта в ее граничной точке х—0. С учетом соотношений (17) и (18), опуская индекс «к» порядкового номера скважины, будем иметь:

г„

Р\ = Р

Г х=0 Г

= Р

(х,у)еГд 1 3,д '

Т = т = т

| х=0 ^з,д -'пп'

-2 I х=0 + 0з х=0 '

0|х=0 = 01 |х=0 + 02 01 | х =0 =Р | х =0 р1 , 0з | х=0 =Рз | х=0 (1 -в) Р ,

О21 х=0 = Р2 | х=0 СнР ( Р|х=0, Т|х=0 ) .

Расчет термогидродинамических характеристик смеси в обсадной колонне скважины на участке Бок от ее забоя х—0 до координаты 2э подвески УЭЦН

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

Решения систем (6)-(8) и (9)-(13) в области Бок скважины и в насосном узле Бэцн должны быть связаны на их границе - в точке входа в насосный узел. Условия, сопрягающие давление, температуру и состав трехфазной смеси с учетом ее нагрева при обтекании ПЭД, при запишем в виде:

Р = Р

14=0 1 вх : О 4=0 = Овх

т = т

14=0 вх ' О2 4=0 = О2вх ,

Оз 4=0 = Озвх ,

Ьэ

Твх = Тп +((д-М)УО'Ср', ыа = \лмё4.

/ '=1 0 Здесь Тп - температура смеси в скважине на приеме установки (на нижнем торце ПЭД); Нэ - мощность, потребляемая насосным узлом (мощность на валу двигателя); параметры потока в точке х—2э_0, найденные из решения уравнений (6)-(8) в области Бок, обозначены нижним индексом «вх».

Таким образом, система уравнений

Р = Р х=2Э+0 вых ' Т =т х=2Э +0 вых'

О = о х=1Э + 0 — вых > о2 = о х=2Э +0 _ ^вых '

Оз 1 х=2Э+0 Озвых

позволяют рассчитать термогидродинамические характеристики водонефтегазового потока в НКТ при 2э<х<Ис.

Постановку (1)-(21) будем называть прямой задачей, в которой расчет осуществляется последовательно в областях Бпп, Бок, Бэцн, Бнкт вплоть до устья скважины. Однако в практике нефтедобычи важным эксплуатационным параметром является давление на устье добывающей скважины при х—ИС (давление Ру,дО определяемое по показанию буферного манометра на фонтанной арматуре скважины), которое должно превышать линейное давление Рпин на начальном участке сборного трубопровода для обеспечения транспортировки продукции. Вообще говоря, величина Ру,д(/) может значительно превышать значение Рпин, если для ее регулировки используется специальный штуцер. Поэтому, задачу (1)-(21) целесообразно решать, полагая заданным давление на устье каждой к-ой добывающей скважины:

Р\х=ИС,к Ру,д,к ,

к = 1, М д

1х=ИС ,к у.д,к' ' А

В этом случае постановку (1)-(22) будем называть обратной задачей, в которой величины Рз д к(/), в(^), являются неизвестными, а ее решение должно быть найдено таким образом, чтобы с заданной точ-

ностью согласовать расчетные и заданные значения устьевых давлений Ру,д,к.

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

Итак, в (М+1)-связной области Dnn с внешней границей Г и внутренними границами Гн,, Гд ,k и в областях D0K, D^H и DHKT требуется при 0<t<Tp найти решение системы взаимосвязанных уравнений (1)-(11) с начальным условием (14), условиями сопряжения (15)-(21) и граничным условием (22). Для решения прямой задачи используется метод конечных разностей, а для обратной - итерационные методы.

Соответствующие численные модели и алгоритмы реализованы в среде Visual Studio 2010 на языке C# с использованием библиотеки параллельных вычислений Task Parallel Library (TPL) [15] из платформы .NET Framework 4.0 и набора библиотек NVIDIA CUDA [16,17] в программном комплексе RSWP-2D-II, позволяющем выполнять расчет взаимосвязанных процессов во всех элементах единого комплекса «пласт - система УЭЦН-скважин» с одновременной визуализацией результатов. Специальный программный блок разработан для имитации работы контроллеров станции управления типа «ЭЛЕКТОН».

Распараллеливание расчетов проводится по трем основным направлениям. Во-первых, при решении фильтрационной задачи в неоднородном пласте методами, разработанными в [1,2,18,19], на сетке, которая может содержать десятки и сотни тысяч узлов. Во-вторых, при расчете термогидродинамических процессов в системе добывающих скважин, оборудованных УЭЦН, на основе моделей (6-8), (9)-(13) с помощью численных методов и алгоритмов, аналогичных изложенным в работах [5,9,11,13]. Эти расчеты требуют большого объема вычислений каждом временном слое, поскольку месторождение может быть вскрыто десятками и сотнями скважин, в каждой из которых число Ыд с k узлов при дискретизации соответствующих областей D0K, D^H и DHKT может составлять порядка нескольких тысяч, общее число узлов

N

_ мд

д , с, общ ^ Nд,с

k=1

во всех скважинах может быть как

сопоставимым с количеством сеточных узлов Мплс в пласте, так и значительно превосходить его. Кроме того, при итерационном решении обратной задачи (1)-(22) объем вычислений во всем комплексе «пласт - система УЭЦН-скважин» возрастает пропорционально количеству итераций. В-третьих, при визуализации двумерных полей расчетных и исходных характеристик фильтрационного процесса на дисплее, когда количество пикселей достигает нескольких миллионов. Эти особенности общей модели (1)-(22) идеально подходят для распараллеливания решения каждой из указанных задач, позво-

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

Демонстрация результатов

В качестве примера рассмотрим результаты расчетов характеристик разработки одного из реальных нефтяных месторождений РФ, полученные на основе решения обратной задачи (1)-(22). Структура нефтяного коллектора, имеющего переменную по простиранию толщину И, неоднородна по абсолютной проницаемости К и пористости т, см. рис. 1. Пласт вскрыт рядной системой добывающих и нагнетательных скважин (Мд=51, Мн=15), которые изображены на рис. 1 соответственно темными и светлыми кружками. В нижней части рисунков приведены шкала значений двумерной функции и соответствующая ей штриховая палитра, используемые при построении 2Б-карт. Нетрудно видеть, что толщина пласта по его простиранию изменяется более чем в 150 раз, а проницаемость - на два порядка. Свойства нефти и воды в пластовых условиях и прочие параметры общей модели (1)-(22) задавались так (см. [3,5,11,13]), чтобы они в среднем соответствовали условиям разработки девонских месторождений при пластовых давлении Рпл=16 МПа, температуре Тпл=40 °С, связанной и предельной во-донасыщенности 5*=0.3 и 5*=0.78.

1060 Скважины

32 26 21

lö 10 б 0,2

2120

Скважины : ф — нягчетятельньге о — сойий?

0,3 0,22 0,15 0,1

б

Рис. 1 - Распределение толщины Н (а) и абсолютной проницаемости К (б) по простиранию пласта

Добывающие наклонно-направленные

скважины имеют вертикальное окончание и различные инклинометрические профили, отклонение которых от вертикали может доходить до 300 м. Глубина скважин примерно одинакова (1580 - 1670 м),

а

значения Ру,д, ф) их устьевых давлений в граничном соотношении (22) составляют 1-3 МПа. Для подъема продукции в добывающих скважинах используются погружные установки, техническое оснащение которых зависит от условий эксплуатации скважины и отличается своими паспортными характеристиками. Так, в состав УЭЦН входят насосные узлы, скомпонованные ступенями различных типов как отечественного (например, ЭЦНД5-30, ЭЦНД5-125, ЭЦНД5.1-200, ЭЦНА(К)5А-250(Тв)), так и зарубежного производства (Б8-30, Б20-80, Б13-50, 400Р12 и пр.). Оптимальная производительность 0*вк таких ступеней при работе на воде варьируется в широком диапазоне от 30 до 300 м3/сут, напор дНв* к - от 4 до

6,7 м. Количество Мс, к ступеней в сборках насосных узлов меняется от 100 до 250 штук, координата 2Эк к подвески УЭЦН в скважинах - от 70 до 450 м. В зависимости от мощности насосных узлов установки оснащены различными электродвигателями (например, ЭДБ20-103В5, ЭДБ125-103В5 и др.) с номинальной мощностью от 20 до 125 кВт. Число Лгплс сеточных узлов в пласте превышает 13000, в отдельной скважине Ыд с к~2000, во всех добывающих скважинах N й ~100000.

д , с , общ

Рис. 2 иллюстрирует карты двумерных распределений давления Р и водонасыщенности в один из моментов времени на этапе разработки месторождения после прорыва воды в добывающие скважины.

lnfin 212П 31ЯГ1 _ 424П

~ — нагнетательияе о — öwiiearviyüe

111) ...........Ill

170 135 120 115 110 105 100 90 60

1060 2120 3180 4240

Скважины :  — хагхетатедтяе О — добшающие

0,7 0,65 0,6 0,55 0,5 0,4 0,3 б

Рис. 2 - Распределение расчетных полей давления Р (а) и водонасыщенности S (б) по простиранию пласта

Сплошными линиями на рис. 2а нанесены изобары в окрестности скважин, а на рис. 2б - гра-

ница между областью законтурной воды и нефтеносной залежью. В силу неоднородной структуры пласта условия работы скважин существенно отличаются: обводненность 0к их продукции варьируется от 10 до 80 процентов, суммарные дебиты Qk - от 30 до 280 м3/сут, величины забойного давления Рз д к -от 60 до 100 МПа.

Как показал анализ результатов расчетов, модель (1)-(22) качественно и количественно верно отражает основные особенности процессов при течении многофазных потоков во всех областях Впл, -0ок, _Оэцн и Днкт, существенно отличающихся как своими размерами, так и строением. Она позволяет не только оценить влияние тех или иных термогидродинамических факторов на характеристики движущихся в пласте, лифтовых трубах и насосном узле смесей, но и обосновать возможность и пути повышения эффективности работы отдельных элементов механизированного подъемника и всей эксплуатационной системы в целом.

Перейдем теперь к обсуждению результатов параллельных и последовательных вычислений, приведенных в табл. 1 и табл. 2.

Таблица 1

№ Intel Core

i5-3337U, 2,6 GHz, 2 ядра i7-860 2,8 GHz, 4 ядра i7-X98 3,6 GH 6 яде 0 z, 3

Tn, сек Tlb сек Kll Tn, сек Tlb сек Kll Tn, сек Tlb сек Kll

1 28,8 23,7 1,21 31,2 25,3 1,23 26,7 20,8 1,28

2 114 61 1,87 107 46 2,32 93 35 2,66

2' 123 70 1,76 114 53 2,15 98 40 2,45

3 485 246 1,97 445 170 2,62 397 124 3,2

4 90 37 2,4 78 18 4,34 69 11 6,28

Таблица 2

Intel Core

i5-3337U i7-860 i7-X980

NVIDIA GeForce

GT 640M 500 MH 384 яд LE, z, ра GTX 260, 576 MHz, 192 ядра GTX 58 772 MH 512 яде , z, р

к/с к/с Си Чъ к/с к/с Cll к/с к/с Cll

1,8 3,4 1,9 1,8 4,3 2,4 2 7,2 3,6

9,2 5,1 7,6 4,2 13,8 6,9

В табл. 1 приведены времена Тц и Тп, требуемые для решения тестовых задач на пятидесяти временных слоях соответственно при последовательных расчетах и с использованием различных центральных процессоров Intel. Здесь кц=Тп/Т|| -показатель эффективности параллельных вычислений. Номера вариантов в первом столбце табл. 1

а

соответствуют: 1 - фильтрационной задаче в пласте, 2 и 2' - прямой задаче расчета процессов в системе «пласт - скважины» (2 и 2' - без визуализации и с визуализацией всех результатов расчетов в пласте, скважине и УЭЦН), 3 - обратной задаче для системы «пласт - скважины», 4 - расчетам процессов только в скважинах в прямой задаче. В табл. 2 представлены скорости и ц|| (кадр/сек) визуализации 2Б-полей в пласте в полноэкранном FullHD-разре-шении при последовательном и параллельном режимах на центральных процессорах Intel тех же моделей (верхние значения в ячейках последней строки таблицы) и на графических процессорах NVIDIA с использованием технологии CUDA (нижние значения ячеек). Параметр QpUn/иц характеризует ускорение визуализации в параллельных режимах.

Нетрудно видеть, что в зависимости от типа процессора время Ty параллельного счета при решении задачи фильтрации сокращается на 20-30%, обратной задачи (1)-(22) - в 2-3,2 раза, при расчете процессов непосредственно в скважинах - в 2-6,3 раза, при визуализации результатов - в 2-3,6 раза на центральных процессорах и в 4-7 раз на графических.

Выводы

1. Дана постановка новой модели, описывающей взаимосвязанные процессы тепломассопе-реноса в неоднородном пласте, вскрытом системой нагнетательных и добывающих скважин, оборудованных УЭЦН. Разработаны соответствующие численная и алгоритмическая модели, реализованные в расчетной программе с использованием технологий параллельных вычислений. Показана работоспособность построенной модели.

2. На основе анализа результатов многовариантных вычислительных экспериментов с помощью программного комплекса RSWP-2D-II показана возможность ускорения параллельных вычислений по сравнению с последовательными. При этом коэффициент эффективности параллельных расчетов зависит от типа ПК.

3. Программный комплекс RSWP-2D-II позволяет:

- определить максимальный дебит из пласта на забое конкретной добывающей скважины с учетом состава продукции, соответствующий оптимальным условиям работы насосного оборудования при заданной величине устьевого давления;

- выполнить экспресс-анализ эффективности текущих эксплуатационных режимов скважин, оборудованной УЭЦН конкретного типа;

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

вала ПЭД для оптимизации параметров работы всей

системы добывающих скважин.

Литература

1. А.Н. Чекалин, Численные решения задач фильтрации в водонефтяных пластах. Изд-во Казан. гос. ун-та, Казань, 1982. 208 с.

2. А.Н. Чекалин, В.М Конюхов, А.В Костерин, Двухфазная многокомпонентная фильтрация в нефтяных пластах сложной структуры. Казан. гос. ун-т, Казань, 2009. 140 с.

3. М.А. Пудовкин, А.Н. Саламатин, В.А. Чугунов, Температурные процессы в действующих скважинах. Изд-во Казан. гос. ун-та, Казань, 1977. 168 с.

4. А.Н. Саламатин, Математические модели дисперсных потоков. Казань: Изд-во Казан. гос. ун-та, 1987. 172 с.

5. В.М. Конюхов, Дисперсные потоки в нефтяных скважинах. Изд-во Казан. гос. ун-та, Казань, 1990. 160 с.

6. Р.И. Нигматулин, Динамика многофазных сред. Ч. I. Наука, Гл. ред. физ.-мат. лит., Москва, 1987. 464 с.

7. Л.Р. Минибаева, А.Г. Мухаметзянова, А.В. Клинов, Модели турбулентности для адекватного описания поля скорости в аппаратах с перемешивающими устройствами. Вестник КТУ, 1, 469-478, (2010).

8. М.Х. Хайруллин, М.Н. Шамсиев, П.Е. Морозов, А.И. Абдуллин, В.Р. Гадильшина, И.Т. Салимьянов, Численное решение прямых и обратных задач тепломассопе-реноса в нефтяных пластах. Вестник КТУ, 24, 125-129, (2013).

9. В.М. Конюхов, И.В. Конюхов, ВАНТ. Сер. Матем. моделирование физ. процессов, 4, 60-69 (2012).

10. П.Д. Ляпков, Нефтяное хозяйство, 5, 38-40 (1979).

11. В.М. Конюхов, И.В. Конюхов, С.В. Краснов, ВАНТ. Сер. Матем. моделирование физ. процессов, 3, 34-44 (2012).

12. П.Д. Ляпков, В кн.: Справочное руководство по проектированию разработки и эксплуатации нефтяных месторождений. Добыча нефти. Недра, Москва, 1983, C. 237 - 293.

13. В.М. Конюхов, И. В. Конюхов, С. В. Краснов, В сб. «Супервычисления и математическое моделирование». Тр. 14 Междунар. семинара (Саров, 1-5 октября 2012). ФГУП «РФЯЦ-ВНИИЭФ», Саров, 2013, С. 353 - 362.

14. Технологический регламент ОАО «САМАРАНЕФТЕ-ГАЗ» по эксплуатации УЭЦН. Версия 4.0. Самара, 2007. 41 с.

15. Г. Шилдт, C# 4.0: полное руководство. Пер. с англ. ООО "И. Д. Вильямс", Москва, 2011. 1056 с.

16. J. Sanders, E. Kandrot, CUDA by example: an introduction to general-purpose GPU programming. Addison-Wesly, Upper Saddle River, NJ, 2010. 313 c.

17. А.В. Боресков, А.А. Харламов, Основы работы с технологией CUDA. ДМК Пресс, Москва, 2010. 232 с.

18. В.М. Конюхов, А.Н. Чекалин, М.Г. Храмченков, Миграция разноплотностных жидкостей в водоносных пластах сложной структуры. Изд-во Казанского мате-мат. общества, Казань, 2005. 160 с.

19. Р.Н. Дияшев, Р.С. Хисамов, В.М. Конюхов, А.Н. Чекалин, Форсированный отбор жидкости из коллекторов с двойной пористостью, насыщенных неньютоновскими нефтями. ФЭН, Казань, 2012. 247 с.

© В. М. Конюхов - д.ф.-м.н., проф. каф. прикладной математики К(П)ФУ, V1adimir.Konyukhov@kpfu.ru; И. В. Конюхов -асп. той же кафедры; А. Н. Чекалин - д.ф.-м.н., ведущий научный сотрудник НИЦ им. Н.Г. Чеботарева Института математики и механики им. Н.И.Лобачевского К(П)ФУ, Anato1ii.Cheka1in@kpfu.ru.