УДК 622.276
ПОВЫШЕНИЕ УРОВНЯ КОНТРОЛЯ И УПРАВЛЕНИЯ СИСТЕМ ППД ПОСРЕДСТВОМ СОЗДАНИЯ УНИВЕРСАЛЬНОЙ МОДЕЛИ
IMPROVEMENT OF THE LEVEL OF CONTROL AND MONITORING OF THE RESERVOIR PRESSURE MAINTENANCE SYSTEM THROUGH CREATION OF
THE UNIVERSAL MODEL
С. И. Грачев, А. В. Стрекалов, А. Т. Хусаинов
S. I. Grachev, A. V. Strekalov, A. T. Khusainov
Тюменский индустриальный университет, г. Тюмень
Ключевые слова: техническая гидросистема (ТГС); гидросистемы продуктивных пластов (ГПП); моделирование; алгоритм; комплексная модель Key words: technical hydrosystem (THS); hydraulic systems of productive formations (HPS); modeling; algorithm; integrated model
Большинство ограничений в известных моделях теории гидравлических цепей (ТГЦ) связаны с фиктивными граничными условиями, ограничениями на вид структуры системы, требованиями к виду функций (замыкающих отношений), отсутствием взаимосвязи между гидравлическими режимами и важными техническими показателями элементов (например, взаимодействие с природными системами, перемерзание участков, аварийные режимы работы насосных агрегатов, изменение состояния обратных клапанов, переход гидравлической энергии в тепловую и т. п.).
В данной статье рассматриваются возможности повышения уровня контроля и управления систем ППД посредством создания универсальной модели, с помощью которой стало бы возможно контролировать систему и предсказывать ее поведение при внесении каких-либо изменений. Использование этой модели, прежде всего, позволило бы увеличить энергосбережение систем и точность соответствия технологическим режимам заводнения за счет комплексной оптимизации.
Рассматривается ТГС с произвольной структурой [1, 2], состоящей из т узлов, из которых t узлов являются транзитивными, п звеньев и с путей возможного перемещения текучей среды между активными узлами. Будем считать, что для каждого звена ie[jib, jie] (звена i принадлежащего узлам jib и jie), где jib и jie — его начальный и конечный узлы, задан закон гидравлического воздействия, связывающий перепад давления др (обусловленный техническими свойствами элемента i)
на концах звена и установившийся расход qi:
Др, = f (q,). (1)
Функции f (q¡) — замыкающие отношения характеризуют взаимосвязь перепада давления от расхода, обусловленную внутренними параметрами звена i. Вид f (q¡), например, зависит от параметров гидравлического сопротивления трубопроводной арматуры, производительности насосов и т. д.
Полный перепад давления на концах звена i будет зависеть от функции fi (qi ) и гидростатического перепада, при условии нахождения ТГС в поле гравитации
Др =Ap'i~Az¡ = f (q,) -Az,, (2)
где — гидростатический перепад давления Azi = pg(zu - ze), где р — плотность текучей среды, g — ускорение свободного падения, Zi и z „ — высоты узлов jib и jie над уровнем моря. Влияние факторов «гидростатического парадокса» во внимание не принимается.
Согласно методу «путевой увязки», для любого потокораспределения должны выполняться следующие условия. Во-первых, в каждом транзитивном (соединенным с более чем одним звеном) узле j должен соблюдаться материальный баланс, отвечающий принципу неразрывности (сплошности) потока текучей среды:
^д, = 0, j = 1,2,.../, (3)
где слева стоит алгебраическая сумма расходов по всем звеньям, имеющим общий (независимо от того, конечный это или начальный) транзитивный узел j. Причем если звено входит в узел, то знак перед qi берется положительным, а если выходит - отрицательным. Для активных узлов уравнения материального баланса не записываются.
Во-вторых, сумма перепадов давления др на концах звеньев, входящих в путь
r, должна быть равна сумме гидростатических перепадов давления на концах звеньев, входящих в этот путь, и перепадов давления между узлом начала пути и узлом конца пути. Ими являются активные узлы, символизирующие накопители текучей среды (НТС), давление в которых задано на момент расчета. Для пути r можно записать
1др/ =1 f (я,) =P - P„ +lAz,, (4)
r r r
где слева стоит алгебраическая сумма перепадов давления (обусловленных техническими свойствами объектов) на концах звеньев, входящих в путь r; справа — разность давлений в активных узлах, образующих путь (Prb — давление в узле начала обхода пути, Pre — давление в узле конца обхода), и сумма гидростатических перепадов давления на концах звеньев, входящих в путь.
Направление «обхода» пути задается выбором одного из пары активных узлов начальным, а другого конечным. То есть как и в звеньях, но уже для цепочки от одного НТС до другого. Поскольку активные узлы отражают элемент НТС, то, согласно первому свойству НТС, значения Prb и Pre, характеризующие стабилизированный потенциал текучей среды, должны быть заданы на текущий момент времени (давления в точках возможного притока/оттока — в реках, озерах, емкостях,
пластах и т. д.). Введем вектор q расходов, вектор перепадов давлений Y , обусловленных внутренними свойствами элементов, вектор полных перепадов давлений Y, вектор P давлений во всех узлах модели и вектор Z гидростатических перепадов давлений на концах всех звеньев:
где / — номер звена; п — количество звеньев в структурной схеме; } — номер узла; т — количество узлов в структурной схеме ТГС.
Введем известную в ТГЦ и теории графов матрицу А соединений т узлов и п звеньев, однозначно описывающую структуру системы, безотносительно к кон-
кретным параметрам звеньев и геометрическим координатам узлов. В данной матрице на пересечении строки ], отвечающей узлу ], и столбца г, соответствующего звену г, помещается элемент
0, если звено г не соединено с узлом ]; -1, если звено г исходит из узла]; +1, если узел г является для звена конечным.
Размерность полной матрицы будет т Хп.
Введем матрицу путей В, где на пересечении столбца г, соответствующего звену г и строки г , соответствующей пути г, помещается элемент
0, если звено г не существует в пути г;
+ 1, если звено г входит в путь г и ориентировано согласно обходу этого пути; -1, если звено г входит в путь г и ориентировано против обхода этого пути.
Автором был впервые разработан алгоритм поиска оптимальной матрицы В. Система уравнений в общем матричном виде
А-ё = 0__. (5)
В ■ У ' = Е + и
Здесь вектор Е составлен из разностей давлений РгЬ- Рге между активными узлами соответствующего пути г, а вектор и из гидростатических перепадов между этими активными узлами.
Из матрицы А исключаются строки, соответствующие активным узлам, что необходимо для обеспечения энергетического баланса гидравлических энергий: разность суммы гидравлических энергий в единицу времени, поступающих в активные узлы извне, и суммы гидравлических энергий в единицу времени, исходящих из активных узлов в рассматриваемой гидросистеме, должна равняться сумме гид-
п п
равлических мощностей N звеньев гидросистемы ^ = = X ^ (%) % .
¡=1 ¡=1
Таким образом, имеем систему нелинейных уравнений (СНАУ) относительно неизвестных установившихся расходов — Q. Единственным недостатком исходной системы уравнений (5) является необходимость поиска системы из с линейно-независимых путей. Решение системы (5) осуществляется численным методом Ньютона при нулевом начальном приближении с коррекцией приращений для определения частных производных замыкающих отношений (1) в конечном виде.
В методе «узловой увязки» потокораспределения предлагается новая форма записи системы уравнений и метод ее решения. Так как давления в активных узлах заданы, а приток и отток гидросистемы определяется в процессе решения, система уравнений записывается относительно неизвестных давлений в транзитивных узлах. Основным отличием предлагаемого метода и его математического описания является использование произвольных замыкающих отношений описанных в табулированном виде или полученных в результате решения уравнений моделей элементов, вместо аппроксимирующих зависимостей, которые хотя и удобны в использовании, но являются весьма неточными.
Основой для системы уравнений является материальный баланс в транзитивных узлах, выраженный через зависимости д! = (Ар) расхода в звене г от перепада давления на его концах. Функция 8(Ар) является обратной функции /(<?), то есть для ее определения в произвольной точке — Ар0 необходимо в общем случае решить нелинейное уравнение /(д)-Ар0 = 0 относительно неизвестного расхода Необходимость в численном решении обусловлена тем, что алгебраически отразить обратную зависимость невозможно хотя бы потому, что зачастую замыкаю-
щие отношения /(д) задаются «кусочно» для различных режимов течения. В описываемом методе подразумевается максимально точное нахождение обратных зависимостей.
Выразив неизвестные расходы в (3) через функции д = £ (Ар) и заменив Ар = Р]Л - Р 1„, получим уравнения для t транзитивных узлов, где в каждом уравнении суммируются 81(Ар1) для звеньев, соединенных (смежных) с транзитивным узлом у.
Е ^ (Рь - р1.) = у = (6)
¡е У
Причем давления в транзитивных узлах являются неизвестными, а давления в активных узлах — константами или функциями от времени, которые рассчитываются в модели ГПП.
Для отражения влияния гидростатических перепадов в (6) необходимо добавить элементы вектора г:
Е Ри - р. + **.) = ° у = и, .
¡е у
После приведения (6) к более удобному для решения виду окончательно получим однородную систему нелинейных алгебраических уравнений (СНАУ) относительно неизвестных давлений в транзитивных узлах:
р1; р2,.., Ру , рт ) = ° ^ (Р1, Р2,.., ру , Рт ) = 0, (7)
^ (Р1, Р2,-, Ру , Рт ) = °
где (р1; р2,.., ру,.., рт) — функция зависимости суммы массовых или объемных
расходов потоков, сходящихся в транзитивном узле у , от давлений в смежных с ним узлах (в том числе и активных).
Для описания условий сжимаемости текучей среды необходимо функционально определить зависимость физических свойств среды, влияющих на распределение потоков от давления. Для этого зададимся функциями р(р) — зависимости плотности от давления и у(р) — зависимости кинематической вязкости от давления. Для воды систем ППД автором получены эмпирические зависимости
р(р,Т) = 1000,26-°,°°9Т 1-837-°-°°°2135р + °,43°6р , (8)
, ^ 1,5778 0,00138 (9)
у( p, т)= 0,1846+е^т-
где р — давление, МПа; Т — температура, 0С; р — плотность, кг/м3; у — кинематическая вязкость, мм2/с. Эмпирические константы берутся с соответствующей размерностью.
Здесь необходимы другие замыкающие отношения: функции Б'( р, , р, , 2, , 2, ),
' ¿¡Ь Зхе ^¡Ь Зхе
связывающие массовый расход М, давления и отметки высот концов звеньев. Подставив данные функции в (6), получим систему из уравнений вида
ЕаЛ(р, ,р. ,2. ,2. ) = 0, (10)
где 2у, — абсолютные отметки положения узлов начала и конца звена ' относительно отсчетной плоскости.
Задавшись вектором абсолютных отметок всех узлов, в матрично-векторном представлении (10) будет
Л-У (Р, V) = 0. (11)
Для нахождения зависимостей у'(р р г- г - ) при формировании модели
каждого звена необходимо численно решить уравнение, связывающее массовый расход в звене I с давлениями на его концах. Разделим звено I на N частей. Будем нумеровать каждый участок звена индексом к, начиная от узла начала (рис. 1).
Рис. 1. Схема расчета для условий сжимаемости текучей среды
На каждом малом участке звена ^ ~ N будем полагать величины плотности и
_1_ N
кинематической вязкости постоянными, то есть независящими от изменения давления. Рассматривая функцию зависимости перепада давления на участке д i звена i от массового расхода в звене Mh плотности р и кинематической вязкости v на этом участке в виде функции fi(Mi, р, v, А1), получим следующее соотношение
N
w СM) = Х f (Mi' Pk v>д1, Дги) = Po - PN = Pjb - Pj„, (12)
k=1
где Mi = qkPk — массовый расход в звене, равный произведению объемного расхода на участке k на плотность; Pk = р(pk) — средняя плотность на участке k; Vk = v(pk )
— средняя кинематическая вязкость на участке k; д^ = i0—pg = Д1^—^ pg
— гидростатический перепад давления на участке k.
— Pk ^ Pk 1
Здесь Pk = k 2 k1 — среднее давление на участке k. Давления к началу следующего участка определяются последовательно, на основании замыкающих отношений (1) для несжимаемой ТС:
Pk = Pk-1 - f ^, Pk ,vk, д1) + .
При решении (12) в момент нахождения входящих в (10-11) Si (Pjb, Pjee, zjb, zje) величины P0, z0 и PN, zN являются константами, однозначно определяющими вид зависимости от M; слева. Таким образом, функцию w i (M i) посредством численного решения можно обратить, то есть решить относительно неизвестных давлений в виде функции M = Si(P P z. z. ).
i Jib Jie Jib Jie
Порядок расчета w{ (Mi) зависит от соотношения направления потока и ориентации звена. Так, при противоположной ориентации звена потоку последовательность расчета Pk следует начинать с узла - конца звена, так как причинно-следственная связь прослеживается согласно направлению потока. Вследствие такой неоднозначности, вид функций S'(P P z . z . ) будет несколько изме-
* Jib Jie Jib Jie
няться при последовательном приближении потокораспределения в ходе численного решения (10-11).
Для систем ППД учет теплового обмена важен из-за необходимости учета изменения свойств текучей среды и выявления возможных фактов перемерзания участков ТГС.
Допустим, что для каждого звена 1 известно распределение температуры окружающей среды по длине звена /, описываемое функциональной зависимостью Н(/). Такие зависимости могут быть представлены в произвольном виде: алгебраически, табулированного множества [Нк, /к], в виде констант, интерполяционной зависимостью эмпирических данных и т. д. Предполагается, что теплопередача между текучей средой в звене и окружающей это звено средой происходит под действием перепада температуры потока и окружающей среды и может быть описана для каждого участка А звена исходя из его морфологических свойств (например, площади поверхности контакта сред — у), свойств материала (например, коэффициент теплопередачи — у), свойств текучей среды и перепада температуры между потоком и окружающей средой — А/ в виде функций О (М, А/, у у, А/).
Также предполагается, что нагрев потока обусловлен переходом части гидравлической энергии потока в тепловую («термогидравлический» переход) вследствие гидравлического сопротивления, а также вследствие кинетического воздействия активных элементов насосов на поток. Для трубопроводной арматуры, термогидравлический переход энергии будет описываться исходя из потерь гидравлической энергии в звене 1 на участке А/ для несжимаемой жидкости в единицу времени, как
а (м ,, А/) _ [ Р, А / ] Р, Вт.
Для участка к длиной А/ звена 1 изменение температуры будет складываться из двух составляющих: роста температуры вследствие гидротермического перехода АТ(а,А/)и роста или падения температуры вследствие передачи тепла между
1 к
потоком и окружающей средой АТ(т,д/)
<к
(г) _ ёХМ,,А/) __I Р ) Р __[_Р
[1 -й(»)] Л^р , А/ ] М [1 -й(»)] /, I М , А/
АТ ,
1к м с м ,с, рС,
где Су — удельная теплоемкость текучей среды, кг. К .
Изменение температуры вследствие теплопередачи с окружающей средой на
участке к длиной А/ звена 1 — АТ/т'А/; _ О[М1'М/,у,У'А^], где А/к _ Тк -Н1(/к) —
разность температур потока и окружающей среды в звене на участке к. Таким образом, для нахождения температуры потока на участке к звена необходимо суммировать все приращения температуры на участках с от 1 до к:
к
Т _ Т +У (АТм' + АТ(т,А'
1ьь ¡с >, причем если полагать постоянство вязкости и плот-
ности на участке, то величина А71(а А/) по звену изменяться не будет. Здесь Т1]Ь — температура в узле (он может быть начальным или конечным для звена), в котором поток входит в звено. Для определения температуры на выходе из звена:
N
Т. _ Т. +У (АТ(а'А/) +ДТ.(т'А/)).
¡„ ¡Ь У ¡к ¡к '
к _1
Для трубопроводов без учета гидротермического перехода в зависимости от температуры потока предыдущего участка к
т,= тм + (т - тм )е
у-у-
м,с„
где ТМ — температура окружающей среды звена , на участке к; у — коэффициент
Вт
теплопередачи, мк; у = Ж-А/-(г — поверхность контакта потока и внешней среды для труб круглого сечения, м2.
Таким образом, после замены ТМ на Н,( а /-к) окончательно получим формулу для определения температуры потока в звене - круглому трубопроводу для несжимаемой жидкости на участке к относительно участка к-1:
а/- ж - у - (
М , с
Тк = Н , (А/ - к) + [Тк1 - Н, (А/ -к) ] е М,с +
^-й(О)] * [ •А/; к
рс.,
N
где П(О) — функция, описывающая долю рассеиваемой части гидравлической энергии, которая не переходит в тепловую в зависимости от скорости потока О .
Для расчета комплексного потокораспределения при неизотермическом течении сжимаемой или несжимаемой жидкости необходимо совместить потоко- и теплораспределение в системе. Для этого описанные выше зависимости для каждого звена , интегрируются в функции АТ = в1(МпТ] ), описывающие перепад
температуры АТ, потока между температурой на входе — Т, о и выходе из звена. Например, для трубопровода этой функцией будет
(М»ТЛ.) = £
А/- Ж -у- ( г
Н, (А /-к) + 1 - Н, (А /-к) ]е Мс +
М£„
1 - И
у
А РкЧ,
(13)
Причем здесь при к = О То = Ти, то есть температура на входе в звено.
Во избежание возникновения ситуации бесконечного АTi при М, = О будем полагать, что при нулевом массовом расходе распределение температуры потока по длине звена будет эквивалентно распределению температуры внешней среды, то есть согласно Н (/ ).
Рассмотрим задачу термораспределения при текущем найденном потокорас-пределении. Допустим, после решения (11) имеем распределение давлений р и массовых расходов М для всех узлов и звеньев модели. Также заданы граничные условия термораспределения: температуры в активных узлах, в которых происходит приток (при данном потокораспределении) в гидросистему. На основании зависимостей (13) для каждого звена , возможно рассчитать распределение температуры во всей системе.
Вводятся функции температуры т. (ТХ,Т2,...,Т.,...,Тт), определяющие на основании (13) всех звеньев зависимость температуры в узле . от температуры во всех узлах.
Т] (ТИе] ,Тк2е] ,
Ткпв] )
£1М ,,к| \Тк +в^к (М ^ Т)]
км._
£ м «
кв.
(14)
1
Ч
к = 1
+
где k е у — индексы узлов инцидентных узлу ], из которых в узел у есть приток;
i е к — индексы звеньев, соединяющих узлы к и узел у; п — количество узлов, инцидентных узлу у, из которых в узел у есть приток. Причем задающими температуру в узлеу считаются узлы, смежные с ним, из которых в узелу имеет место приток.
Для нахождения температур Ту в узлах предлагается следующая СНАУ, решаемая методом простой итерации.
XI\тк (M„k Л)]+xlSl
т (T T T ) = —-
'lV* klel'1 k2el>-"> 1knel> n
X M ek + x
т (T T T ) =kem
bm V klem' L k lem'"''1 knem )
XI Mekl\Tk +di£k (M,.Ek ,Tk )]+ XmSm
X M iek + Xm
(l5)
На рис. 2 показаны зависимости распределения температуры по длине звена от массового расхода при H(l) = const = 5 0C.
Рис. 2. Зависимость распределения температуры по длине трубопровода от массового расхода жидкости:
снизу —вверх 10, 20, 50, 100,500, 1100, 3000 т/сут
Как видно из графиков, учет термогидравлического перехода вносит существенные уточнения модели течения: по мере увеличения М,^ проявление термогидравлического перехода увеличивается, а проявление передачи тепла во внешнюю среду падает. Причем существует такое распределение температур, что охлаждение потока со стороны внешней среды будет полностью компенсироваться нагревом потока вследствие гидравлического трения. Для нашего примера при М] от 500 до 1 000 т/сут.
Для использования комплексной модели системы ППД необходимо объединение граничных условий модели ТГС (МТГС) и модели ГПП. С этой целью для МТГС удобно в наиболее простой схеме рассматривать давление в активных узлах в виде зависимости от времени — Р(/), которая будет обусловливаться моделью ГПП, а для модели ГПП в качестве граничного условия удобно задать зависимостью Q(t) приемистости от времени для каждого звена — скважины.
Необходимо также учитывать динамику проводимостей или замыкающих отношений Дд) для призабойных зон пласта (ПЗП) нагнетательных скважин. Следует полагать, что проницаемость ПЗП также должна описываться зависимостью Щ), которая будет рассчитываться в модели ГПП на каждом шаге времени.
При пересчете комплексного потокораспределения — (КПР) во времени будем иметь динамику термогидравлического состояния ТГС в виде комплексных функций КПР(0, зависящих от внутренних свойств ТГС и от состояния ГПП. В работе [1] детально отражена разработанная автором модель ГПП с учетом условий динамической взаимосвязи с рассмотренными здесь моделями ТГС.
Разработанные автором модели были реализованы в виде программного комплекса Ну^а'8ут [4], который был внедрен на ряде нефтегазодобывающих предприятий Западной Сибири.
Список литературы
1. Стрекалов А. В. Системный анализ и моделирование гидросистем поддержания пластового давления. - Тюмень: ИФ «Слово», 2002. - 324 с.
2. Стрекалов А. В. Математические модели гидравлических систем для управления системами поддержания пластового давления. - Тюмень: Тюменский дом печати, 2007. - 664 с.
3. Меренков А. П., Хасилев В. Я. «Теория гидравлических цепей». - Н.,1985. - 276 с.
4. Стрекалов А. В. Свидетельство о регистрации программы для ЭВМ № 2002611864. Комплекс универсального моделирования технических гидравлических систем поддержания пластового давления (Нуёга^ут). -М., 2002.
Сведения об авторах
Грачев Сергей Иванович, д. т. н., профессор, заведующий кафедрой «Разработка и эксплуатация нефтяных и газовых месторождений», Тюменский индустриальный университет, г. Тюмень, тел. 8(3452)283027, e-mail: grachev-si@mail.ru
Стрекалов Александр Владимирович,
д. т. н., профессор, Тюменский индустриальный университет, г. Тюмень, тел. 89220040077, e-mail: HYdRASYM@mail.ru
Хусаинов Артем Тахирович, к т. н., доцент, Тюменский индустриальный университет, г. Тюмень, тел. 89220070151, e-mail: iQ-tyumen@mail. ru
Information about the authors
Grachev S. I., Doctor of Engineering, professor, head of the chair «Development and operation of oil and gas fields», Industrial University of Tyumen, phone: 8(3452)283027, e-mail: grachevsi@mail.ru
Strekalov A. V., Doctor of Engineering, professor, Industrial University of Tyumen, phone: 89220040077, e-mail: HYdRASYM@mail.ru
Husainov A. T., Candidate of Science in Engineering, Industrial University of Tyumen, phone: 89220070151, e-mail: iQ-tyumen@mail. ru