УДК 622.451:622.414.416
Б.П. Казаков, д-р техн. наук, зав. лаб., aero [email protected] (Россия, Пермь, Горный институт УрО РАН),
A.B. Шалимов, канд. техн. наук, ст. науч. сотр., 8912-48-58-977, [email protected] (Россия, Пермь, Горный институт УрО РАН)
АДАПТАЦИЯ МЕТОДА УЗЛОВЫХ ДАВЛЕНИЙ К РАСЧЁТАМ ВОЗДУХОРАСИРЕДЕЛЕНИЯ В РУДНИЧНЫХ ВЕНТИЛЯЦИОННЫХ СЕТЯХ
Проведён аналитический обзор численных методов моделирования воздухорас-пределения в вентиляционных сетях. Исследованы на сходимость в различных условиях алгоритмы методов последовательных приближении, контурных расходов и узловых давлений. Установлены причины плохой сходимости метода узловых давлений, на основании анализа которых произведена корректировка метода и устранены выявленные в нём недостатки.
Ключевые слова: воздухораспределение, аэродинамические сопротивления, естественная тяга, тепловые депрессии, напорные характеристики, невязка давлений, сходимость.
Современные рудники и шахты представляют собой сложные разветвлённые многоуровневые сети, состоящие из большого количества горизонтальных, вертикальных и наклонных выработок. Проветривание таких вентиляционных сетей осуществляется, как правило, несколькими источниками тяги, - наземными и подземными вентиляторами главного проветривания и вспомогательными вентиляторами местного проветривания. Определённый вклад в движение воздуха вносят также естественная тяга и тепловые депрессии, возникающие в негоризонтальных выработках при значительной разнице температур воздуха в них. Потребность в моделировании проветривания вентиляционных сетей подобной сложности вызывает необходимость в разработке специализированных математических методов расчёта, позволяющих с достаточной точностью и за реальное время определять воздухораспределение в различных условиях.
Вентиляционная сеть задаётся ориентированным графом, аэродинамическими сопротивлениями всех ветвей и напорами всех источников тяги. Математическая постановка задачи расчёта воздухораспределения заключается в составлении уравнений 1-го и 2-го законов сетей. Пусть сеть состоит из n ветвей и v узлов. Первый закон сетей записывается в виде v-1 уравнений в узлах:
о = ’E(±Q'I"), (!)
j
где s - номер узла (от 1 до v-1); j - номера выработок, инцидентных узлу с номером s; Q(jS - расходы воздуха. Поскольку общее число уравнений
должно быть равно числу неизвестных п, то оставшееся число уравнений т=п-(у-1). Это уравнения, отражающие второй закон сетей:
где I - номер (от 1 до т) уравнения (независимого контура); у - номера выработок, содержащихся в контуре с номером I; Аж(/} - напор источника тяги в выработке с номером у контура с номером I.
Схема сети задаётся в виде ориентированного графа, на котором все выработки имеют направление предполагаемого движения воздуха (ветви со стрелкой). Стрелки расставляются без нарушения первого закона сетей во всех узлах, т.е. нет узлов, в которые стрелки бы только входили или бы только выходили. В подсистеме уравнений (1) для определённости можно принять, что знаки перед расходами положительны, если стрелки направлены в узел, и отрицательны, если направлены из узла. В подсистеме уравнений (2) знаки зависят от направления обхода контура: по часовой или против часовой стрелки (выбор направления может быть любой). Если направление напора источника тяги совпадает с направлением обхода,
то ставится знак “+”, если нет - то “-“ (знак же самого напора, так же как и знак расхода определяется совпадением с направлением стрелки выработки). Если направление стрелки выработки совпадает с направлением обхода, то перед ставится знак “+”, если наоборот, то “-“. Знак расхода
в (2) от направления обхода не зависит и определяется только направлением стрелки выработки: “+” - по стрелке, “-“ - против стрелки.
Система уравнений (1) - (2) является нелинейной, и для её решения, вообще говоря, могут быть использованы стандартные математические методы решения нелинейных уравнений, например, метод Ньютона. Однако, как показывает практика решения вентиляционных задач, метод сходится только при достаточно «удачном» выборе начального приближения, подбор которого для систем большой размерности является весьма длительной и трудоёмкой процедурой. Следует заметить, что метод Ньютона предназначен для решения уравнений и систем уравнений произвольной нелинейности, а система (1) - (2) таковой не является. Она состоит из двух подсистем, в которых (1) является линейной, а вторая (2) - квадратичной, т.е. слабо нелинейной. Это обстоятельство послужило поводом к разработке других, более эффективных специализированных методов расчёта систем уравнений такого типа. Первым из них был метод последовательных приближений (МПП), разработанный в первой половине прошлого века [1]. Суть метода заключается в последовательной поконтурной увязке расходов с постепенным уменьшением невязки давлений в контурах. Метод достаточно прост в реализации, изначально был предназначен для ручных вычислений (электронно-вычислительной техники в то время ещё не было)
(2)
]
]
и имеет широкую область сходимости по начальным приближениям в отличие от метода Ньютона. С появлением и развитием ЭВМ метод совершенствовался на предмет ускорения сходимости (путём увязки расходов сразу по нескольким контурам) и на сегодняшний день является наиболее широко используемым при решении сетевых вентиляционных задач.
С наибольшей точностью и математической корректностью процедура уменьшения невязки давлений в контурах была реализована в методе контурных расходов (МКР), в котором увязка расходов производится сразу во всех независимых по расходам ветвях сети [2]. Метод достаточно сложен в реализации, занимает значительно больше оперативной памяти, чем МПП, сходится на порядок быстрее и может по праву считаться наилучшим вариантом метода последовательных приближений. К сожалению, ввиду того, что необходимое для программной реализации метода развитие вычислительной техники было достигнуто не так давно, МКР для решения вентиляционных задач распространения пока не получил.
В работах [2] и [4] описан ещё один специализированный метод решения системы уравнений (1) - (2) того же уровня сложности и математической корректности, что и МКР, но принципиально от него отличающийся. Это - метод межузловых депрессий (ММД) ([4]) (другое название -метод узловых давлений ([2])). Суть метода заключается в уменьшении невязки расходов воздуха в узлах путём увязки падений давлений во всех не -зависимых по давлениям выработках ([2]) (либо увязки давлений во всех узлах ([4])). Таким образом, ММД симметричен МКР относительно порядка решения подсистем (1) и (2). Преимуществом ММД перед МКР является тот факт, что уменьшаемая невязка расходов в узлах позволяет говорить о точности их вычисления, в то время как результат решения уравнений с помощью МКР позволяет оценивать лишь точность вычисления давлений. По результатам проведённых численных экспериментов ([2], [3]) установлено, что ММД имеет гораздо худшую сходимость, чем МКР. В чём же причина такой разницы в сходимости двух аналогичных численных методов? Один метод (МКР) имеет отличную сходимость, другой (ММД) -практически не сходится. Выяснить причину расходимости ММД и найти возможность её устранения можно на основе сравнения и анализа математических алгоритмов обоих методов.
Алгоритм метода контурных расходов имеет следующую структуру.
1. Задаются начальные приближения расходов воздуха Qk в хордах графа (независимых по расходам выработках сети).
2. Из решения линейной подсистемы уравнений (1) находятся значения остальных расходов воздуха в ветвях дерева графа (зависимых по расходам выработках сети).
3. Подстановкой значений расходов воздуха в правую часть подсистемы уравнений (2) находятся невязки давлений во всех независимых контурах.
4. Поскольку падения давлений во всех выработках могут быть представлены как функции независимых расходов, то увязка этих расходов может быть связана с невязкой давлений в независимых контурах с помощью стандартной процедуры разложения функции многих переменных в ряд Тейлора с сохранением слагаемых только первого порядка (по аналогии с методом Ньютона):
І,к
7 , (3)
где I - номер независимого контура от 1 до т;} - номер ветви от 1 до п; к - номер хорды 1 до т; ЬР, - невязка давлений в ¿-м контуре; О- - расход воздуха в }-й ветви; Qk - расход воздуха в к-й хорде; Дл^ Дл^О}) - напорная характеристика источника тяги в }-й ветви (если нет, то 0); Щ - сопротивление }-й ветви в 1-м контуре (-Щ, 0 или в зависимости от наличия ветви в контуре и её ориентации относительно обхода контура). Производ-
О „ (1)
ные —— находятся из решения подсистемы уравнении (1), производные
дОк д( Ап-)
определяются аналитически или численно в зависимости от
»2,
сложности напорных характеристик. Из решения системы уравнений (3) определяются значения увязочных расходов Д2к.
5. Определяется очередное приближение расходов воздуха в хордах 2кг+1 = 21 + АОк (г-номер итерации), и процесс повторяется с пункта
2. Расчёт прекращается, когда значения невязок давлений в контурах становятся меньше заданного значения.
Применительно к ММД система уравнений (1)-(2) записывается в
виде:
0 = £ (+Р<‘>) , (4)
І
1 = Г + ь*7(2,) , (5)
о = Е(Щ”), (6)
где £ - номер узла; I - номер независимого контура;} - индекс суммирования по номерам ветвей, инцидентных узлу с номером б в подсистеме (6), и по номерам ветвей входящих в контур с номером I в подсистеме (4); Р} -
падение давления вместе с напором источника тяги Дл}- в ветви с номером} (с индексом £ - для ветвей, инцидентных узлу s, с индексом I - для ветвей
контура 1). Знаки выбираются также как и в системе (1)-(2). Алгоритм метода узловых давлений имеет следующий вид.
Задаются начальные приближения падений давлений Р£ в ветвях дерева графа (независимых по давлениям выработках сети). Либо задаются давления во всех узлах графа, а Р£ находятся.
Из решения линейной подсистемы уравнений (4) находятся значения остальных падений давления на хордах графа (зависимых по давлениям выработках сети).
По формуле (5) определяются все расходы воздуха.
Подстановкой значений найденных расходов в правую часть подсистему уравнений (6) находятся невязки расходов во всех узлах.
Поскольку расходы воздуха во всех выработках могут быть представлены как функции независимых падений давлений, то увязка этих падений может быть связана с невязкой расходов в узлах по аналогии с МКР:
где £ - номер узла от 1 до V;} - номер ветви от 1 до п'^ - номер ветви дерева 1 до V; §О3 - невязка расходов в ¿-м узле; Р - - падение давления в }-й ветви (с напором источника тяги); Рг - падение давления в ¿-й ветви дерева (с напором источника тяги); - сопротивление }-й ветви (-Щ, да или +Щ в
зависимости от инцидентности ветви ¿'-му узлу и её направленности - в узел (+) или из узла (-)). Если Щ/<0, то знак «-» выносится из под корня.
ар.
Производные —- находятся из решения подсистемы уравнений (4). Из
ЭРг
решения системы уравнений (7) определяются значения увязочных давлений ДРГ.
1. Определяется очередное приближение перепадов давлений в ветвях дерева Р(г+1 ^ = Р[ + АР{ (г -номер итерации), и процесс повторяется с пункта 2. Расчёт прекращается, когда значения невязок расходов в узлах становятся меньше заданного значения.
Алгоритмы МПП, МКР и ММД были реализованы численно и проанализированы на предмет сходимости при расчёте водухораспределения для вентиляционных сетей реальных рудников. Результат - МПП сходится медленно, МКР сходится на порядок быстрее и ММД не сходится вообще. Причина плохой сходимости ММД, как оказалось, заключается вовсе не в преобладании недиагональных элементов в узловой матрице (Максвелла), как это заявлено в [2], а в большой разнице аэродинамических сопротивлений выработок реальных рудников (от 10-10 кмюрг для сбоек до 1000 кмюрг для перемычек). Численный эксперимент показал, что ММД начинает сходиться при разнице сопротивлений не более трёх порядков, а при разнице в два порядка сходится не хуже МКР. Таким образом, причина расходимости
1 дР,
(7)
ММД заключается в том, что при большой разнице сопротивлений двух соседних ветвей на одной из них с большим сопротивлением падает весь напор, а на другой ветви с меньшим сопротивлением падение напора, практически, равно нулю. Ноль оказывается в знаменателе выражения (7), что и не даёт итерационному процессу ММД сойтись, - он осциллирует с огромным шагом вокруг этих «нулевых» решений. Выравнивание же сопротивлений «уничтожает» нулевые решения и позволяет ММД сойтись. Исходя из этого, можно сделать вывод, что метод узловых давлений в классическом виде не пригоден для расчёта воздухораспределения в реальных вентиляционных сетях, поскольку область сходимости метода в три порядка на разброс сопротивлений выработок слишком ограничена.
Несмотря на указанные недостатки, ММД может быть модернизирован на предмет улучшения сходимости следующим образом. Чтобы избежать появления нулевых величин в знаменателе, в качестве независимого базиса следует взять не перепады давлений, а расходы воздуха в ветвях
ДО], а не по АР-. В результате этих преобразований изменяются пункты 5 и 6 алгоритма ММД (пункты 1-4 остаются прежними):
5. Поскольку расходы воздуха во всех выработках могут быть представлены как функции расходов в ветвях дерева, то увязка этих расходов может быть связана с невязкой расходов в узлах:
где ^ - номер узла от 1 до V; і - номера инцидентных узлу б ветвей; 1 - номер ветви дерева 1 до V; Ь0$ - невязка расходов в я-м узле; Оі - расход ві'-й
шения подсистемы уравнений (4), после чего определяются производные
литически или численно в зависимости от сложности напорных характеристик источников тяги А,кі = А%і(Оі) . Из решения системы уравнений (8) определяются значения увязочных расходов ДО,.
6. Определяется очередное приближение перепадов давлений в
дерева, т.е. не Р, а 2і
и раскладывать невязки расходов в узлах по
(8)
ветви; Ог - расход в ,-й ветви дерева. Производные —- находятся из ре-
дР,
ер
определяются ана-
ветвях
дерева Р,(г+' > = Н (0(’г + &в, )\0,г> + Д0,| - Ап, (о(г> + АО,)
(r-номер итерации), и процесс повторяется с пункта 2. Расчёт прекращается, когда значения невязок расходов в узлах становятся меньше заданного значения.
Проведённая корректировка алгоритма себя полностью оправдала -ММД начал сходиться при любых сопротивлениях ветвей (правда, с несколько меньшей скоростью, чем МКР), и может использоваться в расчётах, как альтернатива методу контурных расходов.
В заключение следует отметить другой способ улучшения сходимости ММД, предложенный авторами [4], основанный на введении дополнительного линейного (ламинарного) сопротивления выработки Ял и заменой стандартного квадратичного закона AP = RQ на комбинированный
AP = R^ + RQ . Эффективность такого способа борьбы с расходимостью метода очевидна, поскольку введение дополнительного линейного слагаемого решает «проблему деления на ноль» в области малых расходов. Однако физическая корректность подобного подхода к моделированию аэродинамического сопротивления горной выработки представляется сомнительной, т.к. в реальных условиях даже медленное движение воздуха остаётся турбулентным, а, значит, и второй закон сетей остаётся квадратичным.
Список литературы
1. Меренков А.П., Хасилев В.Я. Теория гидравлических цепей. М.: Наука. 1985. С. 280,
2. Шкундин С.З., Иванников А.Л. Разработка метода межузловых депрессий для расчёта вентиляционных сетей в нормальных и аварийных условиях // Труды научного симпозиума «Неделя горняка -2010». Горный информационно аналитический бюллетень. Вып. 1. 2010. С. 448-458.
B. Kazakov, A. Shalimov
Adaptation of nodal pressure to calculating air-distribution in mine ventilation jets
The state-of-the-art review of numerical methods of modeling distribution of air streams in ventilating networks is lead. Algorithms of methods consecutive approach, plani-metric charges and central pressure are investigated on convergence in various conditions. The reasons of bad convergence of a method of central pressure on the basis of which analysis updating of a method is made are established and the lacks revealed in it are eliminated.
Key words: air-distribution, aerodynamic resistance, natural draft, heat depression, forces characteristics, residual of pressure, convergence.
Получено 22.09.10