CC BY
25
Формирование пространственной структуры поля концентрации в многорукавном русле*
В, А, Шлычков
Институт водных и экологических проблем, СО РАН, Новосибирск, Россия
e-mail: [email protected]
The paper presents a numerical model of 2D currents in channel systems of a complicated spatial configuration. The model is based on the two-dimensional Saint-Venant equations. The influence of inhomogeneous in the flow structure on the admixture transport processes in water was investigated. A calculation of the flow characteristics and passive substance redistribution was presented for the case of the 16-km section of the Ob river near Novosibirsk.
Введение
Исследование закономерностей перераспределения субстанции в речных руслах играет важную роль при анализе и интерпретации различных физических процессов, таких как перенос загрязняющих примесей в воде, формирование ледовых и шуговых полей в период ледохода, перемещение твердых взвесей и русловые деформации и др. Особенности транспорта концентрации в речном потоке обусловлены локальными динамическими условиями, которые в значительной степени связаны с морфометрией руслового ложа.
Речные системы Западной Сибири отличаются сложным многорукавным строением водотоков с многочисленными протоками, островами, поймами и узлами деления потока. Установление количественных связей в многосвязной системе русел представляет собой в этом случае серьезную теоретическую задачу, а компьютерная реализация численных моделей для расчета потоков сложной топологии основана на использовании нетривиальных вычислительных подходов. Как правило, применяются одномерные уравнения Сен-Венана с формулированием условий сопряжения потоков в точках ветвления.
Методы решения одномерных систем на клеточных комплексах (графах) изложены в [1], где изучается общая краевая задача для больших систем дифференциальных уравнений, когда областью определения служит не отдельный отрезок, а совокупность отрезков, образующих граф. Решение задачи о движении воды в системе русел осложняется тем, что при взаимодействии водотоков течение на одном из участков влияет на другой. Поэтому не всегда возможно независимое назначение граничных режимов во внутренних узлах сети. При этом моделирование многосвязных (кольцевых) областей
* Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований (грант № 05-05-98012)
© Институт вычислительных технологий Сибирского отделения Российской академии наук, 2008.
в некоторых случаях чревато появлением ложных решений, приводящих к накоплению фиктивных циркуляций вокруг внутренних петель графа [2].
Значительное внимание при разработке численных моделей в системах водотоков уделялось описанию переноса примесей [3] (акцент при постановке задач делался на правильной формулировке внутренних граничных условий в узловых точках), изучению возможности применения эффективных схем интегрирования, созданию алгоритмов оптимизации. Переносу примесей в многорукавных дельтах рек посвящена работа [4], в которой обсуждаются способы определения коэффициентов диффузии и дисперсии.
Вместе с тем одномерно-графовая схематизация для описания топологически сложных водных объектов имеет недостаточную универсальность, не позволяет учесть эффекты извилистости и меандрирования русел, углов сопряжения водотоков, влияния пойменных массивов, которыми, как правило, перемежаются рукава и протоки. Это обстоятельство приводит к необходимости учета плановой структуры течения. Очевидным преимуществом двумерной модели служит отказ от построения графа узлов ветвления и сложной системы прогонок по ребрам с итерационным согласованием в вершинах.
Цель данной работы заключается в изучении влияния морфологии разветвленного русла на процессы миграции пассивной субстанции в водной среде с применением численной модели плоских (плановых) течений.
1. Постановка задачи
Для расчета гидравлических параметров речного потока используется двумерная модель плановых течений, В горизонтальной плоскости введем декартову систему координат с осями х, у. Поверхность руслового ложа зададим уравнением г = гь(х, у), где гь — функция, описывающая рельеф дна. Уравнения плановых течений запишем в виде [5]
дки дкии дКпу д(к + гь) 9 , ,
-Ж + — + ~а^ = - с||и|“’
дки дкиу дкуу д(к + гь) 9 I I
1Г + “аГ + "аТ = ~9>1~аг~ ~ с!1 и|1''
дк дик дик „
«+аГ+вГ = 0- ()
Здесь £ — время; к — глубина потока; и,ь — компоненты средней по глубине горизонтальной скорости; д — ускорение силы тяжести; С3 — коэффициент Шези; |и| = \/и2 + V2 — модуль скорости течения.
Сформулируем краевые условия. На входном створе, который позиционируем на некотором поперечнике русла, будем считать известным суммарный речной расход В поперечнике выходного створа задается уровень свободной поверхности, пересчитан-
к
условия па компоненты скорости и = V = 0 и пространственное распределение глубин к в момент £ = 0, Искомые параметры течения были получены интегрированием уравнений (1) по времени до выхода на установившийся режим.
Базовые алгебраические соотношения получены на основе метода конечных объемов, Использовались криволинейные сетки с узлами, разнесенными по граням элементарного пространственного бокса [6], Система уравнений (1) проецируется в непрямоугольный элементарный бокс со сторонами, близкими по направлению к криволинейным границам русла. Такой подход обеспечивает детализацию поперечной структуры потока на сравнительно небольшом числе узлов, расположенных на продольных образующих, приблизительно повторяющих изгибы русла. Использовался метод конечных разностей с применением неявных алгоритмов. Пространственная дискретизация дифференциальных операторов основана на представлениях о схемах с невозрастанием полной вариации (Total Variation Diminishing, TVD), гарантирующих монотонность решения путем использования перестраивающегося шаблона и подходящего выбора аппроксимации производных в различных участках численного решения. Монотонизация схемы проводилась по методу Куранта—Изаксона—Риса [7], Неявная часть TVD-операторов алгоритмизировалась в соответствии с методом [8],
2. Результаты расчетов
Сформулированная плановая модель водотоков апробирована на 16-километровом участке нижнего бьефа р, Оби под Новосибирском, Река на данном участке является судоходной, глубины в отдельных местах превышают 11 м. Геометрия русла в плане с элементами криволинейной сетки и распределением глубин показана на рис, 1 (направление течения показано стрелками). Русло имеет сравнительно сложное морфологическое строение: крупный остров Медвежий (обозначен символом О) делит поток на два основных рукава, от которых ответвляются мелкие протоки, образующие группы островов. Ширина рукавов не превышает 800-900 м, по более полноводному левому рукаву проходит судовой ход, В правобережье расположена городская зона отдыха (За-ельцовский парк), ниже по берегам и на островах находятся многочисленные садовые общества.
Цифровую модель рельефа формировали по результатам натурных съемок, проведенных экспедицией МГУ на данном участке реки в 2003-2004 гг, [9], Натурные данные послужили основой для создания цифровой модели рельефа на криволинейной сетке, пригодной для использования в плановой численной модели.
Рис. 1. Сеточная структура расчетной области и контуры русла р. Оби (серая заливка)
Важным этапом вычислительной процедуры является алгоритм разбиения расчетной области на совокупность элементарных объемов (боксов), на которых строятся дискретные соотношения. Сетки произвольной структуры эффективно адаптируются под плановую геометрию расчетной области и особенности морфометрии,
В настоящее время разработаны и широко используются алгебраические, геометрические и дифференциальные методы генерации сеток в областях произвольной формы. Алгебраические методы обладают рядом достоинств, к которым можно отнести простоту реализации, возможность построения координатных линий заданной конфигурации, высокую скорость автоматической обработки узлов и т, и. Наиболее эффективным в данном случае оказался комплексный подход, использующий алгебраический и дифференциально-вариационный методы, обладающий необходимой гибкостью и позволяющий регулировать свойства проектируемой сетки с помощью управляющих функций и параметров. При этом задача формирования системы сеточных узлов в произвольной области ставится для системы квазилинейных уравнений эллиптического типа с граничными условиями, определяемыми соответствием контура области и канонического прямоугольника - прообраза исходной области. Сетка, показанная на рис, 1, имеет 100 х 500 узлов в поперечном и продольном направлениях соответственно и обеспечивает среднее пространственное разрешение 20 х 30 м. Меньший шаг в поперечном к оси русла направлении задан в связи с необходимостью более детального описания изменчивости потока при его перераспределении по рукавам.
Искомые параметры течения получены для заданного расхода воды = 1800 м3/с, соответствующего весенне-летнему периоду, интегрированием динамической задачи по времени до выхода на установившийся режим, В стационарном решении из уравнения неразрывности в (1) следует существование функции тока, которая является характеристикой пространственного распределения расходов. Расчетная конфигурация линий тока показана на рис, 2, Линии тока описывают траектории частиц и вместе с тем образуют трубки тока, в которых значение расхода воды постоянно (на рис, 2 трубки тока
3
ред островом, отражающее неоднородность течения в месте деления потока. Основной водоток проходит по левому рукаву, а правый пропускает не более 20-30 % расхода. Задача моделирования распространения пассивной примеси в двумерной постановке сводится к решению уравнения переноса и диффузии-дисперсии:
дкс дкис д^с д , ^ дс д , ^ дс
~я~Г ^ я I я = я ХЪ ^ ~я~ у7Г~’ ( )
дЬ дх ду дх дх ду ду
Рис. 2. Расчетные линии тока
вде с — безразмерная концентрация, м-3; Ех и Еу — коэффициенты дисперсии. Имея в виду изучение опасных последствий и катастрофических событий на реках, отметим, что на траверзе и. Кудряши (что примерно соответствует положению входного створа расчетной области) расположены городские очистные сооружения, которые осуществляют постоянный выпуск в реку Обь сточных вод. При авариях на станции возможен залповый сброс неочищенных вод высокой концентрации, в связи с чем возникает необходимость количественных оценок попадания загрязнения в парковую часть города по водотоку правого рукава, В соответствии с этим на входной границе заладим постоянное и равномерно распределенное в поперечнике значение с = 1 и рассмотрим процесс формирования поля концентрации по длине русла.
Установившаяся со временем пространственная картина распределения примеси представлена на рис, 3, Можно видеть, что в расчетной области образуется ряд зон
с
ки). Как правило, эти зоны формируются в областях торможения потока, где вынос и рассеяние субстанции затруднены. Отсюда следует, что в случае аварийных сбросов длительного действия можно ожидать аккумуляцию загрязнителей в водотоке правого рукава, заметно превышающую концентрацию источника эмиссии, В правом рукаве с
На рис, 4 показано распределение суммарного по поперечному сечению количества примеси (с) вдоль маршевой координаты русла £. Здесь также хорошо заметна неоднородность концентрации по потоку — значения варьируют в пределах 20-25 % от
Рис. 3. Пространственное распределение концентрации загрязнителя в русле
1000
1200
800
/, кт
600
о
5
10
15
Рис. 4. Пространственное распределение суммарной концентрации загрязнителя
невозмущенного среднего, задаваемого на входном створе (штриховая линия). Из сопоставления кривых на рис, 4 следует, что подмеченная неоднородность концентрации в русле обусловлена неравномерностью распределения скоростей течения: в однородном потоке со скоростью, не меняющейся вдоль русла, концентрация остается постоянной по всей длине, как это показывает штриховая кривая.
Список литературы
[1] Воеводин А.Ф., Шугрин С.М. Методы решения одномерных эволюционных систем. Новосибирск: Наука, 1993. 368 с.
[2] Эвотт М.Б. Гидравлика открытого потока. М.: Энергоатомиздат, 1983. 272 с.
[3] Васильев О.Ф., Воеводин А.Ф. Математическое моделирование качества воды в системах открытых русел // Динамика сплошной среды: Сб. научн. тр. /АН СССР. Сиб. отд-пие. Ин-т гидродинамики. 1975. Вып. 22. С. 73-88.
[4] Мак-Доуэлл Д.м., Коннор Б.А. Гидравлика приливных устьевых рек. М.: Энергоатомиздат, 1983. 312 с.
[5] Стокер Дж.Дж. Волны на воде. М.: Изд-во иностр. лит., 1959. 618 с.
[6] Шлычков В.А. Численное исследование разрешения неустойчивости Кельвина—Гельмгольца в русловом потоке с поймой // Сиб. журн. вычисл. мат. 2007. Т. 10, № 4. С. 417-428.
[7] Куликовский А.Л., Погорелов Н.В., Семенов А.Ю. Математические вопросы численного решения гиперболических систем уравнений. М.: Физматлит, 2001. 606 с.
[8] Harten A. On a class of high resolution total-variation-stable finite-difference schemes // SIAM J. of Numer. Anal. 1984. Vol. 21, N 1. P. 1-23.
[9] Отчет по НИР “Выполнить комплексные исследования гидрологического и руслового режима р. Оби (Новосибирская ГЭС — с. Дубровино), оценить современное состояние и разработать прогноз деформаций русла и посадки русла”. М.: МГУ, 2005. 206 с.
Поступила в редакцию 14 марта 2008 г.
