Математическое моделирование движения водного потока

в нижнем бьефе (трехмерная модель)

М.В. Землянникова – канд. техн. наук, доцент

ФГОУ ВПО «Московский государственный университет природообустройства»,

г. Москва, Россия

В.А. Фартуков – канд. техн. наук, доцент

ЗАО «Бюро сервиса и эксплуатации», г. Москва, Россия

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

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

Модель турбулентного тепломассопереноса основана на уравнении Навье – Стокса (закон сохранения импульса), уравнении неразрывности (закон сохранения массы жидкости), закона сохранения энергии, уравнения диффузионного переноса (закон сохранения массы), к-е модели турбулентности, уравнения переноса «концентрации жидкости в воздухе» для цели аппроксимации свободной поверхности [2, 3, 9].

Моделирование трехмерного течения воды в рассматриваемом створе основано на конечно-объемном методе решения уравнений гидродинамики с применением прямоугольной адаптивной сетки и локальным измельчением [4, 5, 6].

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

Моделирование проводилось для трехмерной турбулентной модели потока воды с применением к–е модели турбулентности, переноса скалярных величин и их флуктуаций. Отображение результатов вычислений осуществляется в виде векторов на плоскости или поверхности течения воды. Интегрирование параметров течения жидкости по сечению и поверхности с отображением в виде изолиний, тоновой заливкой.

Процесс расчета течения воды в нижнем бьефе ГЭС производился в следующей последовательности:

cоздание расчетной области (геометрии объекта), то есть построение цифровой модели нижнего бьефа ГЭС;

задание математической модели;

установка граничных условий;

установка расчетной сетки и критериев адаптации ее (сетки) по граничным условиям;

установление параметров метода расчета, определяющего скорость сходимости алгоритма;

выбор шага по времени расчета;

визуализация результатов расчета, определение и сохранение числовых значений характеристик параметров течения;

оценка точности расчетов методом сходимости по сетке.

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

Входное сечение состоит из трех рассматриваемых участков:

донное водопропускное, здесь задается граничное условие, где скорость втекания воды в расчетную область V ≥ 0;

поверхностный водосброс – граничное условие заключается в задании трех компонент скорости втекания воды в область расчета (V_x, V_y, V_z) – так как в конце конструкции водосброса установлен носок схода под углом 35⁰;

выход из турбин, – здесь граничное условие такое же, как и для донного водовыпуска, то есть значение скорости втекания в расчетную область V ≥ 0.

На выходе из расчетной области. То есть в конце рассматриваемой модели нижнего бьефа ГЭС – ставится условие свободного выхода, иначе - нулевой поток и ∂P/∂n = 0, а давление на границе Р = 0.

Составленная математическая модель движения воды в нижнем бьефе ГЭС представляет собой совокупность уравнений конвективно-диффузионного переноса, которые решаются методом конечных объемов. Производится интегрирование рассчитываемых переменных по объему каждой i-й ячейки расчетной сетки и отрезку времени. Ячейка расчетной сетки имеет форму произвольного многогранника. Решение конвективного – диффузионного переноса осуществляется способом восстановления рассчитываемой переменной в ячейке.

Трехмерное восстановление решается линейной комбинацией трех одноразмерных функций потока жидкости вдоль осей координат X, Y, Z внутри ячейки. Для расчета были применены два способа вычисления:

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

Эта схема расчета была применена на предварительных этапах расчета с целью получения первого приближения решения, оценки общего времени расчета и объема вычислений, с последующим расчетом схемой высокого порядка точности. Последующие расчеты проводились схемой высокого порядка точности ступенчатой функцией, принимающей только два значения min и max во всей расчетной области.

Вычисления проводились явным методом расчета, то есть интегрирование проводилось шагом по времени t, значение которого определяется из условий устойчивости вычислительного алгоритма.

Шаг по времени в явном алгоритме определяется условием Куранта-Фридрихса-Леви

t ≤ t_min = X_ki/V_ki.

где t_min – min шаг, по времени определяемый в результате обхода всех расчетных ячеек; X_ki, V_ki – линейный размер к-й ячейки и скорость в направлении – Х.

Предварительное определение явного шага по времени осуществлялось с помощью числа Куранта-Фридрихса-Леви.

Оценка точности полученного результата вычислений весьма сложна и потребовала большого объема времени по отладке математической модели, которая заключалась в следующем:

проверялась сходимость расчета по сетке, то есть проводилась серия расчетов при одних и тех же начальных условиях на сетке, которая последовательно сгущалась по всей области расчета. Уменьшение расчетных ячеек увеличивало точность решения пропорционально их размерам, с применением порядка аппроксимации расчетной сетки n = 2. Применение этого способа оценки точности расчета позволило оценить точность решения исходных уравнений. Заданная точность вычислений составляла для объема – 0,01, для рассчитываемых характеристик – 0.01.

В основе вычислительного процесса лежит стандартная модель к-е турбулентности, основанная на уравнениях Навье – Стокса. Моделирование турбулентной вязкости осуществляется через выражения величины к-е,

где к – турбулентная энергия, м²/с²; е – скорость диссипации турбулентной энергии, м²/с³; μ – турбулентная вязкость, кг/с/м.

Модель движения двухфазной жидкости включает в себя следующие уравнения:

уравнения Навье – Стокса;

уравнения для энтальпии;

уравнения для концентрации;

уравнения для к и е;

уравнения для переноса функций заполнения жидкостью ячейки.

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

На основе имеющейся информации о рельефе местности, характере течения реки, грунте слагающего ложе русла реки, была составлена цифровая модель участка нижнего бьефа, то есть область расчета.

Исследуемый участок нижнего бьефа охватывал всю ширину (2000 м) и по длине составил – 1500 м.

Построение цифровой модели нижнего бьефа ГЭС проводилась, как указывалось раньше, с помощью САПР. Расчет проводился с помощью адаптивной сетки размером 20 х 15 х 1 с общим числом расчетных ячеек порядка 1700000 и количеством вычислений на расчетном шаге более 680000. Среднее машинное время расчета составила более 40 ч.

Результаты проведенных расчетов показали, что в месте слияния водных потоков водопропускных сооружений, наблюдается сильная турбулентность, с образованием многочисленных завихрений с неравномерным распределением скоростей потока по осям координат X, Y, Z. На рисунках 1, 2, 3, 4 показаны векторы скорости водного потока по направлениям осей координат X, Y, Z. Ближе к донной части их абсолютная величина уменьшается, а на уровне свободной поверхности достигает величины начальных скоростей (V₁ = 7,77…8,0 м/с – скорости потока донного водовыпуска при пропуске расхода и уровнях воды в верхнем бьефе 208,0 и 209,5, соответственно, V₂ = 20,0 м/с – скорость на сходе с носка водосброса и V₃ = 1,3…1,9 м/с скорости потока на выходе из турбин).

Общий характер течения воды в нижнем бьефе характеризуется довольно сильной неравномерностью распределения скорости потока как по ширине русла, так и по длине.

По мере удаления от водопропускных сооружений наблюдается выравнивание поля скоростей, векторные значения по осям координат X, Y, Z сглаживаются, что приводит к затуханию циркуляций потока и снижению вихревых зон. Однако вдоль всего левого берега прослеживается незначительная поперечная циркуляция. Абсолютные значения скорости потока в центре поперечного сечения нижнего бьефа находятся в диапазоне 3,2… 3,4 м/с, по мере приближения к левому берегу скорость несколько снижается, достигая величины в 3,0 м/с.

Рис.1. Компонента скорости Vx при Х = 50 – расстояние от водопропускного сооружения	Рис.2. Компонента скорости Vx при Х = 800 – расстояние от водопропускного сооружения
Рис. 3. Компонента скорости Vy при Х = 500 – расстояние от водопропускного сооружения	Рис. 4. Компонента скорости Vz при Х = 600 – расстояние от водопропускного сооружения

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