Кириллов В.В.

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

Введение

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

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

Определим задачу, которую должна описывать и решать модель ГТС. Пусть заданы начальное распределение давления Pi,0 по всей ГТС (где i=1,..,n v пространственный слой) и граничные условия на подаче и отборе газа в систему по времени. Необходимо определить распределение давления Pi,j и расхода Qi,j по ГТС для каждого временного слоя j, где j=1,...,m. Пример ГТС с элементами приведен на рис.1

В основе данного подхода к построению модели ГТС лежит предположение, что сеть состоит из ряда объектов, взаимодействующих между собой по определенным алгоритмам. Для сети это два класса объектов v узлы и ребра, которые и представляют собственно сеть. На основе элементов этих классов строится сеть необходимой сложности. Основное отличие классов в том, что каждое ребро может быть связано не более чем с двумя узлами, а узел в свою очередь не имеет ограничений по количеству, относящихся к нему ребер.

Класс узлов состоит из следующих типов объектов:

узлы между ребрами (внешние)

узлы по длине трубы, т.е. внутренние узлы ребра.

Исходя из определения, внутренние узлы связаны с двумя соседними частями ребра. Так как вычисления по ним ведутся по математической. модели трубы, примем их как единое целое с объектом ?труба¦.

Класс ребер состоит из следующих типов объектов:

подкачки, т.е. объекты подачи газа в ГТС;

отборы, т.е. объекты отбора газа из ГТС;

трубы;

компрессорные станции (КС).

Рассмотрим каждый тип подробнее.

Подкачки осуществляют подачу газа в ГТС, имеют ссылку только на один узел ГТС (узел, куда осуществляется подача газа). По времени для них задается изменение давления Pi,j или расхода Qi,j (граничные условия), где j=1,-,m - это количество временных слоев для расчетов.

Отборы осуществляют отбор газа из ГТС, имеют ссылку только на один узел ГТС (узел, куда осуществляется отбор газа). По времени для них также задается изменение давления Pi,j или расхода Qi,j (граничные условия), где j=1,...,m - это количество временных слоев для расчетов.

Трубы осуществляют передачу газа по ГТС, имеют ссылку на два внешних узла ГТС (откуда идет поступление газа и откуда идет отбор газа из трубы). Для них задается начальное распределение давления газа Pi,0 в момент времени t0, где i=1,...,n v это количество внутренних узлов трубы.

КС осуществляют увеличение давления между двумя внешними узлами ГТС с сохранением массового расхода газа за счет увеличения потенциальной энергии газа. Другими словами, КС поддерживают по времени определенный массовый расход в определенном направлении. Для полной характеристики КС нам необходимо учитывать физические характеристики КС, затраты газа на внутренние нужды КС и др. Но для упрощения модели КС можно задавать только изменение расхода газа Q0,j по времени, которое по сути будет содержать в себе все остальные параметры КС. Таким образом, КС имеют ссылку на два внешних узла ГТС (откуда идет поступление газа и откуда идет отбор газа из трубы) и для них задается изменение расхода газа Q0,j по времени, где j=1,...,m - это количество временных слоев для расчетов.

Для всех ребер кроме труб изменение давления или расхода по времени задается. Поэтому основная сложность заключается в том, чтобы определить состояния труб и внешних узлов ГТС для каждого временного слоя.

Модель движения газа по трубе

По начальному состоянию газа в трубе и краевым условиям на концах трубы необходимо определить конечное состояние газа на определенном временном слое. Параметрами состояния являются давление газа Pi,j и его расход Qi,j в точке трубы. Все остальные физические параметры газа, трубы принимаются константами на текущий момент расчета. Тогда состояние трубы - это набор точек вдоль трубы (их количество n). Обычно расстояние между ними берут одинаковым (D X). Так как от начального до конечного состояния происходит какой-то промежуток времени T, его тоже можно разбить на m промежутков D t, чтобы можно было проследить изменение состояния в n точках трубы на каждом из m промежутков времени (или временных слоев). В итоге решение - это плоскость распределения давления и расхода по времени и длине трубы.

Приводим возможные варианты состояния газа в трубе:

Стационарное состояние газа (без движения). Когда давление газа одинаково по всей трубе P=const и соответственно расход газа Q=const.

Стационарное движение газа. Когда давление газа неодинаково распределено по всей трубе и рассчитывается по формуле стационарного режима, но расход газа Q=const по всей трубе.

Нестационарное движение газа. Когда давление газа P также неодинаково распределено по трубе, но при этом расход газа Q распределен по трубе также неодинаково.

Исходя из постановки задачи имеется начальное состояние трубы в виде массивов Pi,0 и Qi,0 (где i=1,.,n). Тогда совокупность всех временных слоев (т.е. состояний на текущий момент) описывается состоянием труб Pi,j и Qi,j (где i=1,....,n и j=1,...,m).

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

 (1)

где P - давление газа;

Q - расход газа;

K - коэффициент физических параметров трубы и состояния окружающей среды

Проанализировав данную систему, отметим следующее:

Правило 1. Скорость изменения давления по времени t равна скорости изменения расхода газа по расстоянию x.

Правило 2. По расстоянию x скорость изменения квадрата давлений равна квадрату расхода газа.

Необходимо принять во внимание, что во временном слое можно вычислить расход газа по соседним давлениям. А при перемещении по времени можно определить давление, если известно предыдущее состояние соседних по x расходов.

Построим разностный аналог приведенных выше правил.

На рис.3,а) отображено общее графическое представление производных правила 1. Для построения используется равномерная шахматная разностная сетка по расстоянию x и времени t. Для нестационарного режима течения газа справедливы оба правила 1 и 2 одновременно, а при стационарном режиме правило 1 вырождается, так как изменений по времени для давления и расхода нет и соответственно уравнение правила 2 представляет правило распределения давления и расхода по трубе для стационарного режима вне зависимости от времени.

В идеале каждый временной слой представляется в виде отрезков трубы со стационарным режимом между узлами. Тогда для каждого отрезка соблюдается правило 2. На рис.3,а) представлена схема по правилам 1 и 2. В виде разностной формулы это:

(2)

На рис.3,b) представлена схема расчета давления и расхода газа исходя из правила 1 при переходе с одного временного слоя на другой. Необходимо напомнить, что правило 1 используется только в разрезе одного временного слоя.

На рис.3,а) правило 1 справедливо для точки i,j, но для вычислений такая схема неудобна, так как необходимо знать P и Q на двух временных слоях, чтобы вычислить любой третий. Точно также можно рассуждать и о пространственных слоях. Пространственный слой v это ряд значений по P и Q для точки трубы в разрезе времени. Основная особенность пространственного слоя в том, что он определен либо по P, либо по Q, в то время, как во временном слое, определены P и Q для разных точек. Таким образом, чтобы вычислить значения на одном из трех пространственных слоев необходимо знать значения на двух оставшихся.

На рис.4 представлена схема для расчета пространственного слоя (значения Pi,j+1) с использованием только одного предыдущего временного слоя, а не двух как на рис.3. Для этого шаг на временных и пространственных слоях надо принять одинаковым, что позволит совместить значения и в точке i,j v будет точка их равнозначных значений. Хотя все перечисленное можно применить и для случая, когда и соответственно. Если же, то.

Через точку i,j проходит один временной и один пространственный слой, они пересекаются под углом 90¦ . По правилу 1 для этой точки задействованы окружающие ее 4 точки: (i-1,j),(i+1,j),(i,j-1),(i,j+1). Причем, зная значения трех узлов по временному или пространственному слою, всегда можно найти оставшиеся два узла по соседним слоям. Схемы, приведенные на рис.3,b), построены с использованием этого свойства.

Если мы знаем во временном слое Qi-1,j,Qi+1,j и Pi,j, то

 (3)

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

В идеале, любое нестационарное состояние газа стремится перейти в стационарное. Первоначально v это стационарный режим течения газа, а затем стационарное состояние газа (т.е. выравнивание давления в замкнутом объеме). Тогда в постановке нашей задачи возможны два случая:

Задано распределение давления газа по узлам трубы. Необходимо определить состояние в узлах после определенного промежутка времени без внешних воздействий. Это выравнивание давления в замкнутом пространстве.

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

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

Для формулы (3) v это разница между расходами газа для узла i,j по длине трубы, т.е.

 (4)

где

K- - это приведенный коэффициент физических параметров трубы при изотермическом движении и без учета перепада высот на концах трубы;

Sign(Qi,j,Qi+1,j) - это функция, определяющая знак изменения давления в узле i,j на текущем шаге в зависимости от направления и значения окружающих его расходов газа Qi,j и Qi+1,j. Алгоритм выбора отображен на рис.5.

На рис.5 приведен пример, когда Pi-1,j>Pi+1,j (хотя для и Pi-1,j<Pi+1,j правила описанные ниже сохраняются) причем Pi,j может иметь 5 разных значений (1-5) относительно уровней 1-3. Уровень 1 представляет собой значение Pi,j при стационарном течении газа между Pi-1,j и Pi+1,j, уровни 2,3 v это давления на концах трубы.

По рис.5 видно, что увеличение или уменьшение давления в узле i,j зависит от направления и значения окружающих его расходов газа. Для вариантов 4 и 5 давление в узле i,j экстремально, т.е. в четвертом случае Pi,j>Pi-1,j>Pi+1,j, а для пятого случая Pi,j<Pi-1,j<Pi+1,j. Тогда, если направления расходов Qi-1,j, Qi,j совпадают, то при направлении расходов в узел давление увеличивается (вариант 5), при направлении расходов из узла давление уменьшается (вариант 4). Естественно, если расходы равны Qi-1,j=Qi,j по абсолютному значению и противоположны по направлению (сколько газа поступает столько и уходит), то изменения давления не происходит.

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

Для узла трубы давление увеличится, если в него поступает газа больше, чем отбирается. Давление уменьшится , если газа поступает меньше, чем отбирается.

Теперь определим общую схему алгоритма расчета каждого временного слоя всей сети.

Алгоритм расчета состояния временного слоя ГТС

Сначала просчитывают все трубы ГТС. Порядок значения не имеет. Для них определяется давление в узле следующего временного слоя по формуле (3). Необходимо отметить, что в крайних узлах трубы, т.е. на ее концах, давление не считается.

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

Как видно из рис.6,а) для определения давления следующего временного слоя трубы используется два соседних расхода газа около этого узла. Если исходить из этого, то для узла ГТС нужно выбрать какое-то ребро с расходом Qi, а все остальные расходы сложить с учетом знака в Qi+1. Тогда давление в узле ГТС можно считать по тем же правилам (3), что и для узла трубы. Тем самым мы обобщаем правило расчета давления для любого узла системы:

Давления узла считается по формуле (5):

 (5)

где - это сумма расходов газа для всех ребер, относящихся к  этому узлу с учетом знака (направления течения газа);

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

Граничные условия для ГТС

Как следует из постановки задачи выше, граничными условиями являются либо задание распределения расхода газа по времени на КС Q(t), либо задание P(t) или Q(t) на подкачках и отборах. Так как для подкачек и отборов задается либо давление, либо расход, то оставшееся неизвестное необходимо определять на каждом временном слое.

Представим, что на одном конце трубы задано P(t), а на другом Q(t). Если мы знаем на каждом временном слое Q(t) для подкачки или отбора, то давление вычисляется по обычной схеме (рис.3). Для P(t) дело обстоит несколько сложнее. На рис.7 рассмотрены оба таких примера.

Из рис.7 следует, что неизвестный расход определяется по известному давлению на следующем временном слое, давлению и расходу текущего слоя.

Заключение

Для доказательства сходимости процесса вычислений по приведенной выше математической модели ГТС были проведены следующие эксперименты. Для трубы с фиксированными физическими параметрами и заданным начальным состоянием по давлению (стационарное течение газа между концами трубы ) наблюдалось выравнивание давления в замкнутом пространстве, то есть при условии, что оба конца трубы закрыты. Каждый эксперимент заключался в изменении числа n (кол-во узлов вдоль трубы на текущий промежуток времени) и числа m (кол-во узлов по времени в определенной точке трубы). В результате была получена сетка m x n значений давления в точке плоскости пространства и времени. На рис.8 значения m и n изменялись в диапазоне 3,....,60. Разброс значений давления P для исследуемой точки (3000сек. на 50м) в полученной сетке максимально равен 1300 Н/м2  (абсолютная погрешность), что составляет около 2% (относительная погрешность) от среднего значения давления в этой точке.

Исследование данной сетки привело к следующим выводам:

Увеличение числа m слабо влияет на увеличение точности вычислений при n>30. В этом случае относительная погрешность составила менее 1% Хотя при 3<n<30 относительная погрешность может достигнуть до 2%.

Увеличение числа n значительно сказывается на относительной погрешности только при m>13 (рис.8b). Однако рекомендуется использовать n>30, т.к. при этом давление в точке равномерно без скачков стремится к своему предельному значению (рис.8d). Природа колебания значения давления около определенной кривой объясняется апроксимацией значения давления разностной схемой решения на грубых сетках. Поэтому лучше использовать n>30, но при этом значительно увеличивать n не имеет смысла., так как увеличится время расчетов.

Необходимо отметить изменение количества вычислений в зависимости от n и m, то есть трудоемкость вычислительного процесса. На рис.9(a-d) представлена поверхность изменения количества вычислений. На ее основе и используя результаты исследований (рис.8) можно сделать вывод об оптимальном сочетании качества и количества вычислений для трубы. Наиболее оптимальными являются значения n=30, m=30 или около этого диапазона. Хотя, конечно, для детализации процесса по времени можно увеличивать значение m, но существенного улучшения относительной погрешности это не даст. Ее изменение составит не более 1%.

http://www.laboratory.ru/articl/math/ram01a.htm

http://www.laboratory.ru/articl/math/ram01b.htm

Список литературы

Сергованцев В.Т., Кучин Б.Л., Гарляускас А.И., и др. Централизованный контроль и оптимальное управление на магистральных газопроводах. vЛ.: Недра, 1973.

Карманов В.Г. Математическое программирование. vМ.: Наука, 1975.

Фильчиков П.Ф. Справочник по высшей математике. vКиев: Наукова думка, 1974.

Мантуров О.В. Курс высшей математики. vМ.: Высшая школа, 1991

Дьяконов В.П. Справочник по алгоритмам и программам на языке Бейсик для персональных ЭВМ. vМ.: Наука, 1987.

Новоселов В.Ф., Гольянов А.И., Муфтахов Е.М. Типовые расчеты при проектировании и эксплуатации газопроводов. vМ.: Недра, 1982.

Бобровский С.А., Щербаков С.Г., Гусейн-Заде М.А. Движение газа в газопроводах с путевым отбором. vМ.: Наука, 1972.

Березина И.В. Ретинский В.С. Оперативное управление системами газоснабжения.

Круг А.Г., Сосоулин Ю.А., Фатуев В.А. Планирование эксперимента в задачах идентификации и экстраполяции. vМ.: Наука, 1977.

Боглаев Ю.П. Вычислительная математика и программирование. vМ.: Высшая школа, 1990.