CC BY
93
DOI: 10.5862/JPM.218.3 УДК: [519.6+004.942]:621.643
Е.А. Михайловский, Н.Н. Новицкий
Институт систем энергетики им. Л.А. Мелентьева СО РАН
модифицированный метод узловых давлений
для расчета потокораспределения в гидравлических цепях
при нетрадиционных замыкающих соотношениях
В статье предлагается методика численного решения задачи потокорас-пределения в гидравлических цепях с сосредоточенными параметрами при произвольных замыкающих соотношениях. Представлены основные типы соотношений для законов изотермического установившегося течения среды по отдельным элементам гидравлической цепи, включая традиционные, неявно заданные по расходу и зависимые от давления рабочей среды. К нетрадиционным соотношениям введены формальные условия применимости, обеспечивающие единственность решения в задачах потокораспределения. Для расчета потокораспределения предложен новый модифицированный метод узловых давлений, имеющий большую универсальность в отношении вида замыкающих соотношений, по сравнению с традиционным, и требующий меньших вычислительных затрат, по сравнению с известной методикой двойных циклов итераций. Приведены анализ и алгоритмизация нового метода, численный пример газотранспортной сети и результаты расчета.
гидравлическая цепь, трубопроводная система, задача потокораспределения, математическое и компьютерное моделирование, модель потокораспределения.
Введение
Трубопроводные системы энергетики тепло-, водо-, нефте-, газо-снабжения и т. п. представлены весьма широким спектром развивающихся во времени и пространстве объектов, которые различны по назначению, масштабам, принципам построения и условиям функционирования. Эффективное решение задач управления их развитием и функционированием имеет важное социально-экономическое значение и непосредственно связано с уровнем применения современных методов математического и компьютерного моделирования.
В Институте систем энергетики им. Л.А. Мелентьева СО РАН разработан в рамках развиваемого здесь научного направления — теории гидравлических цепей [1—3] уникальный арсенал методов математического моделирования, расчета и оптимизации, которые потенциально применимы для трубопроводных систем различного типа и назначения. Вместе с тем, при создании
конечных программно-вычислительных комплексов в разных организациях наблюдается дублирование работ по реализации этих методов. Это связано с необходимостью адаптации методов к прикладной специфике и сфере применения. Неизбежно затрачивается много времени и сил на отладку, развитие и сопровождение таких программно-вычислительных комплексов.
Для преодоления этих проблем мы предлагаем два взаимосвязанных пути:
1) переход на концепцию объектно-ориентированного моделирования трубопроводных систем [4, 5], что даст возможность отделить программные компоненты, реализующие общие методы расчета, от отвечающих за специфику трубопроводных систем или решаемой задачи;
2) развитие методов теории гидравлических цепей с ориентацией на эту концепцию.
Это позволит многократно применять однажды реализованные общие методы теории гидравлических цепей в разных
программно-вычислительных комплексах, развивать расчетные компоненты без их перепрограммирования и, в конечном итоге, — повысить оперативность внедрения методов теории гидравлических цепей для разных типов трубопроводных систем, классов решаемых задач и сфер применения.
В рамках изложенной выше концепции в данной статье рассматриваются вопросы унификации методов расчета потокора-спределения в гидравлических цепях (ГЦ) с сосредоточенными параметрами для обеспечения независимости методов расчета от вида соотношений для законов течения среды по отдельным элементам трубопроводных систем.
В связи с этим излагаются следующие результаты исследований:
структуризация основных типов замыкающих соотношений (законов течения) для элементов ГЦ;
алгоритмизация нового метода расчета потокораспределения в ГЦ при произвольных замыкающих соотношениях;
численный анализ эффективности указанного нового метода.
Традиционная модель и задача потокораспределения
Традиционная модель установившегося изотермического потокораспределения в ГЦ включает законы Кирхгофа и замыкающие соотношения
Ax = Q, ATP = y, y = f(x), (1)
где A — [(m -1) x n] -матрица инциденций узлов и ветвей расчетной схемы с элементами = 1(-1), если узел j является начальным (конечным) для ветви i, и a}i = 0, если ветвь i не инцидентна узлу j; x, y — n-мерные векторы расходов и перепадов давления на ветвях расчетной схемы; f (x) — n-мерная вектор-функция с элементами f(xi), i = 1, n, отражающими законы падения давления от расхода (законы течения) на ветвях ГЦ; Q — (m -1) -мерный вектор узловых расходов с элементами Qj > 0 для притока в узел j, < 0 для отбора в узле j и Qj = 0, если узел j — простая точка соединения ветвей; P — (m -1) -мерный вектор узловых давлений.
Задача состоит в определении векторов х, у, Р при заданных матрице А, векторе 0, известном виде /(х,) для , = 1, п и заданном давлении в одном из узлов (Рт), которое для простоты изложения полагается равным нулю.
Известны многочисленные методы и алгоритмы решения данной задачи, однако, как показано в монографии [1], базовыми являются два метода: контурных расходов и узловых давлений. Оба основаны на методе Ньютона, но с предварительным понижением порядка линеаризованных систем уравнений (1).
Так, классический метод узловых давлений (МД) [1] предполагает поиск решения (1) в пространстве узловых давлений и сводится к организации процесса рк+1 = рк +Арк,
на каждой к-й итерации которого поправка аРк отыскивается из решения системы
А/')"1 АТАРк = -ик,
где ик =Ахк-0; хк = у(ук); ^^ ^ АТРк; /' — диагональная матрица частных производных д// дх„ 1=1,п, в точке хк; у — вектор-функция, обратная к / с элементами У, (у,).
Отсюда видно, что вычислительная схема МД не зависит от вида функций / (х,), к которым предъявляются только требования монотонного возрастания для обеспечения единственности решения задачи потокорас пределения [1]. Однако особенности реализации МД связаны со спецификой вида /,(х,) и соответственно д/ / дх; и у.(у,).
В традиционном случае
/ (х,. ) = | х | - у, где 5,. — гидравлическое сопротивление ветви; у > 0 — приращение давления на активных ветвях (например, с насосными станциями) и у = 0 на пассивных (например, для трубопроводных участков).
Здесь мы имеем явные выражения
д/ / сЦ = 25,. | х,. |,
У,(У,) = ^ У, + у | /5 • у, + у).
Ниже рассматриваются два основных случая нетрадиционных замыкающих со-
отношении, начиная со случая неявных функций у. (у,)..
Замыкающие соотношения, неявные по расходу
Для иллюстрации причин разнообразия замыкающих соотношений рассмотрим известную формулу Дарси — Вейсбаха [6] для потери напора в трубопроводе:
I V2
к = X(V) у—, й 2%
где й, I — диаметр и длина трубопровода; % — ускорение свободного падения; V = V (х) — скорость течения жидкости, х — массовый расход.
Здесь коэффициент гидравлического сопротивления X зависит от числа Рей-нольдса Яе^) = Vd / V, где V — коэффициент кинематической вязкости жидкости, который для изотермического течения полагается константой.
Для вычисления X, в свою очередь, имеется множество формул, отражающих назначение трубопровода, его тип, материал внутреннего покрытия, срок службы, режим течения среды (ламинарный, переходный, турбулентный) и т. п. В качестве
примера в табл. 1 приведены некоторые распространенные в России и за рубежом формулы для коэффициента X и его производной по скорости.
Потеря напора на местных сопротивлениях определяется по формуле Вейсбаха [6]:
к = V2, 2%
где в случае регулирующих элементов для определения коэффициента местного сопротивления £ могут потребоваться аппроксимирующие зависимости. Например, для задвижки приемлема следующая зависимость:
С= «(1 - г)"С - Ъ,
где г — степень прикрытия задвижки; а, Ъ, с — коэффициенты, полученные путем аппроксимации (в том числе «кусочной») таблично заданной характеристики.
Обобщая случаи наличия местных сопротивлений (включая регулируемые) на трубопроводе произвольного типа, для ,-й пассивной ветви ГЦ можно записать, что
(X,) = 5,-(X,) | XI | х,,
(2)
где
Таблица 1
Формулы и производные для коэффициента гидравлического сопротивления
Наименование формулы X
Шифринсона [7] о-п (й Г 0
Прандтля — Никурадзе [8, 9] ' й 1-2 1,14 + 21ё й V К ) 0
Кольбрука— Уайта [6] -21Е ( ке , 2,51 [ 3,7 й Re(V )у[Х)\ 10,04 X V [Сзк-1+ ^ 11п(10) + + ]
Альтшуля [6, 8] / , \0,25 0,1114 + 68 1 ^ й Яе^)) " ^ к- +| 68Г ^ й | v | й) й | v | v
Шевелева [10] А. Г Ао + С1т 1 ^ й У^й ) А1тС (Ах + с )1-т й^1+т
Обозначения: ке — коэффициент эквивалентной шероховатости; А0, А1, С, т — коэффициенты, зависящие от материала трубы и режима течения.
5,. (х,) = X,. (х) + £
I
81
рп2^
сопротивление ветви, так как
V; = 4 Ы /(пй,р), у = 9ф,;
р — плотность транспортируемой - — сумма коэффициентов мест-
здесь среды; С-ных сопротивлений. При этом
/х',= = (25 +5х,х,. )| х,-ах.
где
¿5.
5 . = --
х ,1 1
ах.
¿5-
Ж
8/.
¿5 1 ¿V 1 а х , а х; ¿х а V'
4
рп2^
¿х.
па2р
Таким образом, получение /'х { для трубопровода произвольного типа и назначения требует лишь конкретизации XV, (см. табл. 1).
Гидравлические характеристики активных ветвей, в общем случае, задаются алгебраическими полиномами вида
К,
/(х,) = 2 ,х, Г г-1) х,
q=0
когда
П ,=Е ( ^ 1) Gq ,
q =0
х
|(Ь q,,-2)
Однако они могут различаться показателями степеней ЬдП числом к, и составом слагаемых (при принудительном назначении аС[, = 0 для некоторых q). Также может применяться «кусочная» аппроксимация одной характеристики различными полиномами.
Ранее случаи типа (2) относились к классу ГЦ с переменными параметрами, для которых при расчете потокораспреде-ления применима методика двойных циклов итераций [1]. Согласно этой методике, функции 5,.(х,.) пересчитываются на внешнем цикле, а на внутреннем применяются базовые методы узловых давлений или контурных расходов при фиксированных значениях 5 .
Изложенная выше техника работы с производными неявных по расходу зависимостей у = /() позволяет избежать внешнего итерационного цикла, что резко сокращает вычислительные затраты [4].
Вычисление расхода х к по заданному <-» к
перепаду давлений у к в случае неявных функций типа (2) может быть выполнено итерационно, например, с помощью методов простой итерации, где
х =
У.,
'5(х/ )
з1ЕП( у,),
или Ньютона, согласно которому
х.
= x, +
где к — индекс итерации расчета потокорас-пределения, t — индекс итерации отыскания хк, и] = ук - /, (х]) — невязка.
При t = 0 принимается, что х] := х,-1.
Как известно, высокая скорость сходимости метода Ньютона имеет место при наличии «удачного» начального приближения. Сходимость метода простой итерации в данном случае, наоборот, выше вдали от решения. Поэтому оказался эффективным комбинированный метод, состоящий в однократном применении метода простой итерации перед применением метода Ньютона.
Рис. 1 иллюстрирует процессы сходимости трех предложенных вариантов на примере трубопровода при ук = 5,3085 м вод.ст., d = 300 мм, I = 1000 м; коэффициент X определен по формуле Шевелева для новых стальных труб (А1 = 0,0159; А0 = 1; С = 0,684; т = 0,226), при начальном приближении х 0 = 800 кг/с и решении хк = 71 кг/с. Видно, что применение комбинированного метода дает ощутимый выигрыш, а само решение с приемлемой точностью (и] < 10-4 м вод.ст.) получено всего за три итерации.
Таким образом, в случае неявных функций у (у.) вычислительная схема МД остается без изменений, за исключением применения специальной техники на этапах расчета /' и хк.
Рис. 1. Иллюстрация процесса сходимости разных методов: 1 — простой итерации, 2 — Ньютона, 3 — комбинированного
Замыкающие соотношения, зависимые от давления
Данный случай возникает при невозможности сведения гидравлических характеристик элементов трубопроводных систем к традиционному виду у = / (х). Например, паспортная характеристика газоперекачивающего агрегата (ГПА) связывает не разность, а отношение давлений (степень сжатия) е = рР / рь(рР, рь — величины давления до и после ГПА) с объемным расходом. Наиболее простой аналитический вид этой характеристики [11] —
ф(Рр, Рь, х) = $0Рг \Рг\ - Рь |Рь | - р2х|Х = 0, (3)
где — коэффициенты, полученные путем аппроксимации е2(х) паспортной характеристики ГПА; х — производительность, отнесенная к условиям входа.
также в литературе встречаются другие типы зависимостей е(х) или е2(х) для ГПА, полученные при разных степенях аппроксимирующих полиномов. Например, в работе [12] предложена зависимость
ф(Рр, рь, х) =
(
2 А
4р
Рр|Рр| - Рь|Рь| -
2
(
- Р2
х - Рр 2Р2 р
х РР 2р2 р
(4)
= 0.
В процессе исследования также рассматривались зависимости
ф(рр,рь,х) = — -а - — х-
Рр
Рр
а
,х |х| = 0,
Рр \Рр\
ф( Рр, Рь, х) = р0 Рр |Рр I - Рь\Рь\ -- в1 |рр| х - Р2х|х| = 0,
(5)
(6)
где — коэффициенты аппроксимирующих полиномов для е( х).
С целью обобщения рассмотренных случаев гидравлических характеристик (включая традиционные и неявные по расходу) введем в рассмотрение замыкающее соотношение следующего функционального вида:
ф(рр, Рь, х) = о,
(7)
где рр, Рь — я-мерные векторы давлений, соответственно, в начале и конце ветвей ГЦ. Действительно, так как у = р¥ - Рь, то вместо у = / (х), имеем
Рр - Рь - /(х) = ф(Рр , Рь, х) = 0.
На рис. 2, а показаны примеры классических замыкающих соотношений
у = жх |х|, у = sx |х| - У
Рис. 2. Графики традиционных замыкающих соотношений для пассивной (1) и активной (2) ветвей в координатах х, у (а); х, рг (Ь); х, рь (с)
в координатах х, у. Как отмечалось, для гарантированного получения единственного решения задачи потокораспределения необходимо соблюдение требования монотонности замыкающих соотношений [1], которое можно распространить на рассматриваемый случай (7), представив классические соотношения в координатах х, рР и х, ^^ (рис.2, Ь и с).
Как видно из рис. 2, а, требованию монотонного возрастания у(х) эквивалентно одновременное соблюдение требований монотонного возрастания рг() (рис. 2, Ь) и монотонного убывания по рь(х) (рис. 2, с), что эквивалентно соблюдению условий
дх а дх -> 0, -< 0
дрР дрь
на всей области определения функций (7).
На рис. 3 приведены графики зависимостей (3) — (6) для ГПА. Уравнение (5), в отличие от остальных, не удовлетворяет требованию монотонности. Уравнение (3) является более простым, но менее точным. Уравнение (4) имеет более сложную конструкцию для обеспечения монотонности, однако, и для уравнения (6) монотонность можно обеспечить путем введения фиктивной точки в начале характеристики.
Таким образом, уравнение (6) представ-
Рис. 3. Зависимости давлений на входе (а) и выходе (Ь) газоперекачивающего аппарата от расхода жидкости, выраженные различными соотношениями: 1 - формулой (3); 2 - (4); 3 - (5); 4 - (6)
Таблица 2
Аналитический вид функций для расчета х по заданным р, рь
Используемое соотношение Выражение для x(pF, pL)
(3) JH- sign(a), где a = воPfPf- PlIPlI . в2
(4) а . (а Л , 4 — - sign — + b, Ы / 2 Л где н = Ро + 4Т P'lP'l - PlIPLI , A=Jr P?. V 4в2) JP2
(6) (-b -J\D\ - sign(D)) / 2a, где a = -p 2, b = -Pj \pF\, с = p0 pF\pF\- Pl\Pl\, D = b2 - 4 ac.
ляется наиболее приемлемым как по сложности, так и по точности.
Отметим, что для нахождения хк по заданным давлениям в соотношениях (3) — (6) применимы конечные формулы, приведенные в табл. 2.
Модифицированный метод узловых давлений
Традиционная вычислительная схема метода контурных расходов обеспечивает более быструю сходимость, чем метод узловых давлений. Однако в случае (7) возникают трудности исключения узловых давлений при сведении системы уравнений к контурным расходам. Для газопроводных систем применяют специальные методы [13, 14], например, предполагающие решение систем линеаризованных уравнений относительно поправок к неизвестным давлениям и контурным расходам одновременно, которые из-за этого уступают по эффективности стандартному методу контурных расходов [12].
Здесь излагается новый модифицированный метод узловых давлений (ММД), который представляет собой обобщение МД на случай соотношений (7). Как и МД, он обеспечивает поиск решения в пространстве только узловых давлений и базируется на следующей модификации исходной модели потокораспределения (1):
Ax = Q, pF = ATP, pz = - ALP,
ф(PF ' PL *) = 0,
(8)
где Ар, Аь — [(т -1) х п]-матрицы инциден-ций, фиксирующие отдельно начальные и конечные узлы ветвей так, что АР + Аь = А; ф(Рр, Рь, х) — вектор-функция с элементами Ф, (рР,,■, Рь,I, Х1), ¿ = 1)П отражающими произвольные законы течения, включая зависимые от давления.
Пусть на к-й итерации имеется опреде-
пк
ленное значение Р , которому можно поставить в соответствие значения ркр = Ар.Рк, рЬ = -АЬГРк, а также значение хк, удовлетворяющее уравнению
фФ , рЬ, хк) = 0.
Линеаризация уравнений (8) дает
АДхк = -«к, (9)
Дрр = АГГДРк, дрь = -агдрк, (10)
ф'р дркр + ф'ь ДрЬ + ф' Дхк = 0. (11)
где
« = Ахк - б; фРр = ; фрь = ^; фх =
дрР дрь дх
— диагональные матрицы частных производных порядка п.
Выразим Дхк из уравнения (11):
- Дхк = -(ф' )-1(фРр Др^ +фРь ДрЬ).
С учетом соотношений (10) получим:
= -(Ф х УЧфРЛ - ФРХ )АРк. (12)
Подстановка формулы (12) в равенство (9) дает следующее выражение:
А(фх)-1(фрх - ФРА)АРк = ик. (13)
По правилам дифференцирования неявных функций,
(хР¥ ), =
ох.
(
(дрг ),
дФ,-дх,.
( х'рг, ), =
дх,
(фЫ
дФ, чдх, у
V
Л-1
дФ,-
(дрР), ' дФ,- ТТ,
(др1) ,
, г = 1, п.
Отсюда
хРг = -(фх )-1 фрт ; х'п = -(Ф' )-1 ф рь
- диагональные матрицы соответствующих частных производных порядка п.
Тогда систему (13) можно представить в виде
м ьРк = -ик,
где М = Ах'ррАТ - АхРА -ванная матрица Максвелла. Обозначим
модифициро-
П , = (х') Р К , = -(х'ь) п ,= 1, п.
Тогда диагональные элементы М определяются по соотношению
к = Е п, + 2 к,
,е I- ,е
где 1](1++) - множество ветвей, исходящих из узла (заходящих в узел) ].
Внедиагональные элементы ММ определятся следующим образом:
если узлы ] и , соединены единственной ветвью I, то ММ= -к;, когда узел , является конечным для ветви ¡; Мп = -п,-, когда узел , является начальным для ветви ¡;
если узлы ] и , не соединены, то м, = 0;
если узлы ] и , соединены параллельными ветвями, образующими множество I,, то
к = -2 п,- 2к,.
1ЕI]] , е I],
При соблюдении условий
(х^), >0, (х'ь), <0, I = 1, п,
вытекающих из сформулированных ранее требований к монотонности функций
Ф/ ( рТ ,,, Рь,1 , х,), / = 1,п, имеем неравенства
П, > 0, к, > 0, , = 1, п.
Отсюда Мр > 0 при ] = , и М, < 0 при ] Ф , для всех ],, = 1, т -1. При этом
М ^ 2 к м = 1,т -
,=1,, *]
где строгое неравенство имеет место, когда узел ] (или ,) инцидентен ветви, имеющей другим концевым узлом узел т с заданным давлением.
Покажем также, что введенная модель (8), а также рассмотренная модификация МД накрывают канонический случай замыкающего соотношения у = / (х), когда
Ф(рт, РЬ, х) = рР - рь - /(х) = 0.
В этом случае
Фрт = Еп, Фрь = -Еп, где Еп - единичная матрица порядка п.
При этом
(Фртат - Фра) = (ЕпАТ + ЕпАТ) = = Атр + Аьт = АТ.
С учетом того, что ф = - /х , вме ст о вы -ражения (13) имеем равенство
А/)-1 Ат АРк = -ик,
что совпадает с основной формулой для расчета направления шага в классическом МД [1].
Алгоритмизация метода
Приведем сопоставление вычислительных схем ММД и МД (табл. 3).
В приведенных алгоритмах подчеркнуты те составляющие выражений, которые зависят от специфики замыкающих соотношений.
Таким образом, при объектно-ориентированном подходе компоненты, которые реализуют общий метод решения задачи потокораспределения, обращаются к компонентам, отвечающим за модели элементов. Это обращение дает возможность
Таблица 3
Этапы одной итерации работы алгоритмов метода узловых давлений и его модифицированного аналога
№ Шаг Формула по алгоритму
базовому модифицированному
1 Задано некоторое приближение к решению по давлениям рк - -
2 Вычисление вектора расходов (хк) из формулы АтРк - /(хк) = 0 фф ) = 0
3 Вычисление вектора невязок ик = Ахк - Q; если Ци-к|| < ст, то расчет закончен.
4 Вычисление направления шага А( /х')-1 Ат ДРк = -ик [АхРР Ат Ахрр Ат ]ДР — и
5 Вычисление нового приближения рк+1 = рк + д рк Рк+1 = Рк + 8к Др>к
Обозначение: а — заданная точность соблюдения небалансов расходов в узлах.
получать значения производных и расходов на ветвях с конкретными законами течения в точке текущего приближения.
Компоненты диагональных матриц х'рь, которые необходимы для расчета направления шага по ММД, могут быть по-
Таблица 4
Примеры функции фи ее производные
Вариант функции ф Эф дРF дф дРт Эф Эх
Классический тип замыкающих соотношений
р¥ - рт - зх х + У 1 -1 -2з х
Тип, неявно заданный по расходу
рР - рт - з(х)х х 1 -1 -(2 з + ^х х) |х|
\Pf\Pf - \Рт Рт - з(х)х х 2| РF| -2| Рт| -(2з + з'хх) |х|
Тип, зависимый от давления (на примере ГПА)
Формула (3) 2Ро Ы -2| Р F | -2р2 х
Формула (4) 2| РF | -01 'р в + 4Р2. х - р1 ^ 2Р2 / 21Р F | -2р2 х - Р1 ^ 2 2р2
Формула (6) 2Ро|РF - 21Р F | -РЦ - 2р2|х|
лучены на основе аналитических выражений для
Эф,- дф,- дФ,-Эх,- ' (дРЕ ), ' (Фт ), ' Примеры таких выражений для разных
типов замыкающих соотношений (традиционных, неявных по расходу и зависимых от давления) представлены в табл. 4, где индекс ветви опущен.
Классический МД основан на итерационной формуле
Рк1 = Рк +АРк,
которая не предполагает регулировку длины шага. При этом процесс сходимости часто имеет пилообразную форму и не всегда позволяет найти решение за приемлемое число итераций.
В связи с этим были проведены исследования методов регулировки длины шага 8к в соотношении
Р к+1 = Рк +5кдРк.
В случае строго выпуклой функции невязки || ик(8к) ||2 оказалось эффективным применение длины шага, определенной как
8* =|
и
к-1
(|| и
к-1 112
+1| ик (1) ||2)-1,
где || и ||2 - квадрат евклидовой нормы вектора невязок.
В общем случае более надежным для поиска длины шага оказалось применение метода золотого сечения.
Численные исследования
Для апробации ММД взят пример газотранспортной системы [12] (рис. 4), состоящей из девяти узлов, четырех ГПА и шести пассивных ветвей (трубопроводов). Все ГПА описаны соотношением (4), трубопроводы - соотношением
\рр\рр - |рь|рь = 5х |х|.
Отборы в узлах (млн. м3/сут): ( = 19,1; ( = 14,8; ( = 0,632; ( = 0,32. Фиксированное давление Р9 = 33,778 ат.
Коэффициенты si для трубопроводов: 51 = 0,349; 52= 1,332; 54= 0,436; 55= 4,757; 56= 0,253; 510= 0,006.
Обозначим в, = {р0,, р1;,, р2,}. Тогда в данном примере для ГПА:
в3= {1,040975262; 0,4520492230;
0,1660378943};
в7= {1,049124727; 0,3668417249;
7 0,1867004063};
в8= {1,056105913; 0,4352722015;
0,2396158372};
в9 = {1,040975262; 0,4520492230;
0,1660378943}.
Начальные приближения по давлениям генерировались датчиком случайных равномерно распределенных чисел в диапазоне 0 - 100. Выполнялось 100 таких генераций, и для каждой из них применялся предлагаемый алгоритм. Точность решения проверялась по условиям
||ик|| < 0,01 и ||дРк|| < 0,01.
На рис. 5 представлены результаты такого тестирования в форме процента расчетов, в которых задача была решена за конкретное число итераций . Результатом вычислений являются расходы
х = {19,1; 21,6; 10,8; 14,8; 12,93; 13,25; 13,25; 13,25; 10,8; 2,5}
и давления
Рис. 4. Схема газотранспортной системы:
узлы - окружности с номерами, ветви - линии со стрелками, ГПА - окружности со стрелками на ветвях
Рис. 5. Количество расчетных условий, при которых задача решена за к итераций без регулировки шага (а) и с его регулировкой(Ь)
Р = {31,55; 33,51; 41,76; 32,05; 33,50; 43,80;
44,30; 38,77}.
Проведенные численные исследования на широком диапазоне начальных приближений показали, что ММД сходится во всех случаях. Вычисление оптимальной длины шага позволяет резко (в 4 — 5 раз) сократить среднее число итераций (до 7 ± 1) (см. рис. 5, Ь), что сопоставимо по показателям со сходимостью традиционных методов расчета потокораспределения при классических замыкающих соотношениях.
Заключение
Сформулирована проблема обеспечения общности методов теории гидравлических цепей в сфере реальных информационно-вычислительных технологий для моделирования установившихся изотермических режимов трубопроводных и гидравлических систем различного типа и назначения. Предложено два взаимосвязанных пути ее решения:
1) переход на концепцию объектно-ориентированного моделирования трубопроводных систем, которая обеспечивает
независимость расчетных методов от специфики законов течения среды по отдельным элементам системы;
2) разработка новых методов расчета потокораспределения при произвольных законах течения.
Выделены основные типы соотношений для законов изотермического установившегося течения среды по отдельным элементам гидравлических цепей, включая традиционные, неявные по расходу, зависимые от давления. Для нетрадиционных соотношений введены формальные условия применимости, которые обеспечивают единственность решения задачи потокорас-пределения.
Предложен новый модифицированный метод узловых давлений, который обладает универсальностью в отношении вида замыкающих соотношений, требует меньших вычислительных затрат, по сравнению с известной методикой двойных циклов итераций, может выступать в роли базового метода в рамках концепции объектно-ориентированного моделирования гидравлических цепей.
СПИСОК ЛИТЕРАТУРЫ
[1] Меренков А.П., Хасилев В.Я. Теория гидравлических цепей. М: Наука, 1985. 278 с.
[2] Меренков А.П., Сеннова Е.В., Сумароков С.В. и др. Математическое моделирование и оптимизация систем тепло-, водо-, нефте- и газоснабжения. Новосибирск: ВО «Наука», Сибирская издательская фирма, 1992. 407 с.
[3] Новицкий Н.Н., Сеннова Е.В., Сухарев
М.Г. и др. Гидравлические цепи. Развитие теории и приложения. Новосибирск: Наука, Сибирская издательская фирма РАН, 2000. 273 с.
[4] Новицкий Н.Н., Михайловский Е.А. Объектно-ориентированное моделирование гидравлических цепей // Вестник ИРГТУ. 2012. № 7. С. 170-176.
[5] Рамбо Дж., Блаха М. ИМЬ 2.0. Объектно-
ориентированное моделирование и разработка. 2-е изд. СПб.: Питер, 2007. 544 с.
[6] Альтшуль А.Д. Гидравлические сопротивления. 2-е изд., перераб. и доп. М.: Недра, 1982. 224 с.
[7] Шифринсон Б.Л. Основной расчет тепловых сетей. Теория и методы расчета. М.; Л.: Госэнергоиздат, 1940. 188 с.
[8] СНиП 41-02-2003. Тепловые сети. Госстрой России. Введ. с 1 сентября 2001 г.
[9] Манюк В.и., Каплинский Я.и., Хиж Э.Б. и др. Наладка и эксплуатация водяных тепловых сетей. 3-е изд. М.: Стройиздат, 1988. С. 432.
[10] Актуализированная редакция СНиП 2.04.02-84*. Водоснабжение. Наружные сети и сооружения. Издание официальное. М.: Изд. Министерства регионального развития Российской Федерации, 2012.
[11] Немудров А.Г., Черникин В.и. Расчет режимов работы газопроводов методом определения оптимальных характеристик турбонагне-
тателей // Газовая промышленность. 1966. № 3. С. 31-34.
[12] Евдокимов А.Г., тевяшев А.Д. Оперативное управление потокораспределением в инженерных сетях. Харьков: Вища школа, 1980. 144 с.
[13] Левин А.А., таиров Э.А., Чистяков В.Ф.
Расчет потокораспределения в энергоустановках как гидравлических цепях с регулируемыми параметрами // Трубопроводные системы энергетики: математическое моделирование и оптимизация; под. ред. Н.Н. Новицкого. Новосибирск: Наука, 2010. С. 115-124.
[14] Корельштейн Л.Б., Пашенкова Е.С. Опыт использования методов глобального градиента и декомпозиции при расчете установившегося неизотермического течения жидкостей и газов в трубопроводах // Трубопроводные системы энергетики: математическое моделирование и оптимизация; под. ред. Н.Н. Новицкого. Новосибирск: Наука, 2010. С. 103-115.
СВЕДЕНИЯ ОБ АВТОРАХ
МиХАИЛОВСКМй Егор Анатольевич - младший научный сотрудник лаборатории трубопроводных и гидравлических систем Института систем энергетики им. Л.А. Мелентьева Сибирского отделения Российской академии наук.
664033, г. Иркутск, ул. Лермонтова, 130 [email protected]
НОВицКиИ Николай Николаевич - доктор технических наук, главный научный сотрудник Института систем энергетики им. Л.А. Мелентьева СО РАН, зав. лабораторией трубопроводных и гидравлических систем.
664033, г. Иркутск, ул. Лермонтова, 130 [email protected]
Mikhailovsky E.A, Novitsky N.N. A MODIFIED NODAL PRESSURE METHOD FOR CALCULATION OF THE FLOW DISTRIBUTION IN HYDRAULIC CIRCUITS FOR THE CASE OF UNCONVENTIONAL CLOSING RELATIONS.
In the paper, we propose a method for numerical solving the problem of the flow distribution in hydraulic circuits with lumped parameters for the case of random closing relations. The conventional and unconventional types of relations for the laws of isothermal steady motion of fluid in the individual hydraulic circuit components are considered. The unconventional relations are represented by those implicitly given by the flow rate and dependent on the pressure of the working fluid. In addition to the unconventional relations, the formal conditions of applicability were introduced. These conditions provide a unique solution to the flow distribution problem. A new modified nodal pressure method is suggested. The method is more universal in terms of the closing relation form as compared to the unmodified one, and has a lower computational cost as compared to the known technique of double iteration cycles. The paper presents an analysis of the new method and its algorithm, gives a calculated example of the gas transportation network, and its results.
hydraulic circuit, pipeline system, flow distribution, mathematical and computer simulation, flow distribution model.
REFERENCES
[1] A.P. Merenkov, V.Ya. Khasilev, Teoriya gidravlicheskikh tsepej [Theory of hydraulic circuits]. Nauka, Moscow, 1985.
[2] A.P. Merenkov, E.V. Sennova, S.V. Sumarokov, et al., Matematicheskoe modelirovanie i optimizatsiya sistem teplo-, vodo-, nefte- i gazosnabzheniya [Mathematical modeling and optimization of heat, water, oil and gas supply systems]. VO «Nauka», Sibirskaya izdatel'skaya firma, Novosibirsk, 1992.
[3] N.N. Novitsky, E.V. Sennova, M.G. Sukharev, et al., Gidravlicheskie tsepi. Razvitie teorii i prilozheniya [Hydraulic circuit. Development of the theory and applications]. Nauka, Sibirskaya izdatel'skaya firma RAN, Novosibirsk, 2000.
[4] N.N. Novitsky, E.A. Mikhailovsky, Ob'ektno-orientirovannoe modelirovanie gidravlicheskikh tsepej [Object-oriented modeling of hydraulic circuits], Vestnik IRGTU. 7 (2012) 170-176.
[5] J. Rumbaugh, M. Blakha, UML 2.0. Ob'ektno-orientirovannoe modelirovanie i razrabotka [UML 2.0. Object-oriented modeling and development]. 2nd ed. Piter, St. Petersburg, 2007.
[6] A.D. Al'tshul', Gidravlicheskie soprotivleniya [Hydraulic Resistances]. 2nd ed. Nedra, Moscow, 1982. 224 s.
[7] B.L. Shifrinson, Osnovnoj raschet teplovykh setej. Teoriya i metody rascheta [The basic calculation of heat networks. Theory and methods of calculation]. Gosenergoizdat, Moscow, Leningrad, 1940.
[8] SNiP 41-02-2003. Teplovye seti [The thermal networks]. Gosstroj Rossii. Since September 1, 2001.
[9] V.I. Manyuk, Ya.I. Kaplinskij, E.B. Khizh, et al., Naladka i ekspluatatsiya vodyanykh teplovykh setej [Commissioning and operation of the water heating networks]. 3rd ed. Strojizdat, Moscow, 1988, p. 432.
[10] Aktualizirovannaya redaktsiya SNiP
2.04.02-84*. Vodosnabzhenie. Naruzhnye seti i sooruzheniya. Izdanie ofitsial'noe [The updated edition of SNiP 2.04.02-84*. Water supply. External networks and facilities. Official publication]. Izd. Ministerstva regional'nogo razvitiya Rossijskoj Federatsii, Moscow, 2012.
[11] A.G. Nemudrov, V.I. Chernikin, Raschet rezhimov raboty gazoprovodov metodom opredeleniya optimal'nykh kharakteristik turbonagnetatelej [Calculation of the operating modes of pipelines through the method of determining optimum performance turbochargers], Gazovaya promyshlennost'. 3 (1966) 31—34.
[12] A.G. Evdokimov, A.D. Tevyashev, Operativnoe upravlenie potokoraspredeleniem v inzhenernykh setyakh [Operational management flow in engineering networks]. Vishcha shkola, Kharkov, 1980.
[13] A.A. Levin, E.A. Tairov, V.F. Chistyakov, Raschet potokoraspredeleniya v energoustanovkakh kak gidravlicheskikh tsepyakh s reguliruemymi parametrami, Truboprovodnye sistemy energetiki: matematicheskoe modelirovanie i optimizatsiya [Calculation of the flow distribution in power plants as hydraulic circuits with adjustable parameters, Energy piping systems: mathematical modeling and optimization], ed by N.N. Novitsky, Nauka, Novosibirsk, 2010, 115-124.
[14] L.B. Korel'shtejn, E.S. Pashenkova, Opyt ispol'zovaniya metodov global'nogo gradienta i dekompozitsii pri raschete ustanovivshegosya neizotermicheskogo techeniya zhidkostej i gazov v truboprovodakh, Truboprovodnye sistemy energetiki: matematicheskoe modelirovanie i optimizatsiya [Experience in the use of the global gradient methods and decomposition in the calculation of steady non-isothermal flow of liquids and gases in pipes, Energy piping systems: mathematical modeling and optimization], ed. by. N.N. Novitsky, Nauka, Novosibirsk, 2010, pp. 103-115.
THE AuTHORS
mikhailovsky Egor A.
Melentiev Energy Systems Institute of Siberian Branch of the Russian Academy of Sciences
130, Lermontov St., Irkutsk, 664033, Russian Federation
novitsky Nikolai N.
Melentiev Energy Systems Institute of Siberian Branch of the Russian Academy of Sciences
130, Lermontov St., Irkutsk, 664033, Russian Federation
© Санкт-Петербургский политехнический университет Петра Великого, 2015
