CC BY
30
Вычислительные технологии
Том 12, Специальный выпуск 4, 2007
ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ
ТУРБУЛЕНТНЫХ ТЕЧЕНИЙ МЕТОДОМ КРУПНЫХ ВИХРЕЙ
У. С. Абдибеков, Б. Т. Жумагулов, Б. Д. Сурапбергенов Казахский национальный университет имени аль-Фараби,
Алматы, Казахстан e-mail: uali@kazsu.kz
This research addresses a development of a numerical algorithm for modeling of turbulent flows using large eddy simulation approach. The authors have developed a method for solving the 3D Poisson equation for pressure field using a combination of the Fourier method for one coordinate and the Thomas algorithm for the 2D Fourier coefficients. Several problems have been solved using the proposed algorithm. The results of numerical simulations were compared with the experimental data.
Введение
В реальной природе все течения являются турбулентными. Все они неустойчивые, произвольные и хаотичные в пространстве и во времени. Поэтому исследование турбулентных течений является актуальной проблемой. К наиболее известным методам численного моделирования турбулентных течений относятся прямое численное моделирование (DNS), решение осредненных уравнений Навье — Стокса (RANS) и метод крупных вихрей (LES).
Метод крупных вихрей основан на двух предположениях. Первое состоит в возможности разделения поля скорости на движение крупных и мелких вихрей, причем движение крупных вихрей может быть рассчитано отдельно, что связано с достаточной изотропностью и универсальностью мелких масштабов турбулентного движения. Второе предположение — в возможности аппроксимации нелинейных взаимодействий между крупными и мелкими вихрями только по крупным вихрям с использованием подсеточных моделей. Таким образом, крупномасштабные структуры решаются явно, а мелкомасштабные — моделируются на основе разных моделей и типов функций фильтров [1].
© Институт вычислительных технологий Сибирского отделения Российской академии наук, 2007.
1. Основные уравнения
Основными уравнениями для компонент скорости и являются уравнения Навье Стокса и уравнение неразрывности [2]:
дщ д ^ = др 1 д2щ ¿=12 3- (1)
дЬ дxj г 3 дх^ Р1е дxjдxj'
Для любой переменной / (х) те крупномасштабная часть /(х) получается с использованием функции фильтра С(х):
/(х) = У С (х - х')/ (х') ¿х', (3)
а мелкомасштабная часть /' (х) определяется как
/' (х) = / (х) - / (х). (4)
Существуют различные типы фильтров: • фильтр типа "ящик"
с (|г - г/|) Г1/А, |х - х'|< А/2, С (|Г Г |) \ 0, |х - х'| > А/2;
фурье-фильтр
3
С(|г - г'|) = П
вт^ - х^)/А
(х* - хЭ/А
С(|г - г'|) = (6/пА2)1/2 ехр (-6(х - х')2/А2) .
После применения фильтра к уравнениям (1) и (2) мы получим основные уравнения для крупномасштабных движений:
+ А = + 1 94* _ (5)
дЬ дxj г 3 дх^ Ее дxjдxj дxj'
5 (6)
Ту =щЩ - щЩ, (7)
где тц — иодсеточный тензор, отвечающий за мелкомасштабные структуры, который необходимо моделировать.
Для моделирования подсеточного тензора использовалась модель Смагоринского, которая представляется в виде
3
где ут = С^А^ДцБц)1/2 — турбулентная вязкость; С^ — коэффициент, зависящий от характера течения; А = (А^А^А^)1^3 — ширина сеточного фильтра; Бц — - '
Тц —^-Ткк = —2 ^тБц, (8)
—— ) — величина тензора скоростей деформации |2| дх,- /
2 V дхз
2. Численный метод
Для решения задачи турбулентного течения используется схема расщепления по физическим параметрам [3]:
u* - un
-= - (uraV) u* - z/Au* - Vr; (9)
т
Vu*
Apn+1 = (10)
т
un+l— u*
--— = -Vpn+1. (11)
т
Предлагается следующая физическая интерпретация приведенной схемы расщепления, На первом этапе предполагается, что перенос количества движения осуществляется только за счет конвекции и диффузии. Промежуточное поле скорости находится с использованием метода дробных шагов. На втором этапе по найденному промежуточному полю скорости находится поле давления. Уравнение Пуассона для поля давления решается методом Фурье и представляется в следующем конечно-разностном варианте:
Pi+l,j,k ~ 2Pi,j,k + Pi—l,j,k Pi,j+l,k ~ 2pi,j,k + Pi,j-l,k Ax* Щ
Pi,j,k+1 — 2Pi,j,k + Pi,j,k-1 _ ^ /10ч
+-щ---(12)
Поле давления представляется в виде ряда Фурье [4]:
2 N nkl
Pi,j,k =-fi-L.pl Oi, j,l cos—,
N31=0 N3
/г 2 ^ h nkl Fitjtk = ttZ, Pi hi,icos
N31=0 N3
где
N3 nkl
ai,j,i = Y1 Pk Pi,j,к cos——; k=0 N3
h J? nkl
bi,j,i = L Pk^i,j,к cos——;
k=0 N3
_ ( 1, 1 < i < N3 - 1, Pi 0.5, i _ 0, N3.
После подстановки (13) в уравнение (12) и нескольких арифметических преобразованиях уравнение (12) приводится к следующему виду:
-Ajaj-1 + Bjaj - Cjaj+i _ bj, (15)
где матрицы ai;j;k и bi;j;k представлены в векторном виде aj и bj. Уравнение (15) решается матричной прогонкой с граничными условиями Неймана, После нахождения коэффициентов ai)j,k значения поля давления находятся из формул (13), Для вычисления сумм (13) и (14) применяется .метод быстрого преобразования Фурье, который позволяет вычислить данные суммы за O(N3 ln N3) действий, что сокращает вычислительное время.
На третьем этапе осуществляется ускорение потока под действием сил давления.
(13)
(14)
3. Результаты моделирования
В работе рассматривается модельная задача о каверне. Результаты моделирования сравнивались с экспериментальными работами Prasad и Koseff и с результатами работы Zang, Street и Koseff [5]. Область размерами 1 : 1 : 0.5 (дли на х высотах ширина) покрывалась сеткой размером 50 точек в каждом направлении. Число Рейнольдса, основанное на скорости крышки каверны UB = 1 м/с и длине каверны L = 10.1 м, равно 10 000. На всех границах каверны принималось условие прилипания для вектора скорости. Систе-
L
каверны UB.
Профили средних скоростей U/Uq = u/Ub и W/Uq = w/Uв представлены на рис. 1, среднеквадратичные отклонения пульсаций Urms = ^{и"2)/Ub и Wrms = \J(w"2)/Ub — на рис. 2, напряжения Рейнольдса по центру области UW = {и"w")/U^ — на рис. 3.
Другая задача — это турбулентное перемешивание с вдувом в замкнутой области. Рассматривается экспериментальная область, в которой нагревается нижняя поверх-
Рис. 1. Профили средних скоростей и/Vв Рис- 2. Среднеквадратичные отклонения и уо/Пв пульсаций д/{и"2)/IIв и д/(г///2)Д/в
Рис. 3. Напряжения Рейнольдса по центру Рис. 4. Схема вычислительной области
области {и"п]" )/и'2
0.0 0.2 0.4 0.6 0.8 X
Рис. 5. Значение вектора скорости Ж по середине области (У = 0.502 м)
z
0.8
0.6
0.4
0.2
0.0
41.5 0.3 4.1 0.1 0.3 X
Рис. 6. Значение вектора скорости и по середине области (X = 0.502 м)
б
Рис. 7. Значение поля температуры по середине области: а — Y = 0.502 б — X = 0.502 м
ность, а остальные поверхности имеют температуру вдуваемого воздуха. Размеры области: высота H = 1.04 м, дайна L = 1.04 м и ширина D = 0.7 м. Схема вычислительной области с граничными условиями показана на рис. 4. Высота вдува hin составляет 0.0018 м, скорость вдуваемого воздуха 0.57 м/с и температура вдува и боковых стенок области равна 15 °С. Высота выходного отверстия hout = 0.024 м, температура нагретой поверхности 35 °С. Результаты моделирования сравнивались с экспериментальными данными работы [6] (на рисунках показаны квадратиками). На рис. 5 и 6 приведены сравнительные результаты значений вектора скорости по середине областей W(U3) (У = 0.502 м) и U(U) (X = 0.502 м). На рис. 7 представлены сравнительные результаты значений поля температуры по середине областей.
Заключение
Представленные сравнительные результаты показывают, что построенный численный алгоритм адекватно описывает турбулентные течения в замкнутой области. Некоторые расхождения с численными результатами являются следствием того, что в данной
работе использовались простая модель Смагорипского и однородная сетка во всех направлениях, При решении трехмерного уравнения Пуассона для давления разработан и реализован гибридный алгоритм метода Фурье с применением матричной прогонки для коэффициентов Фурье, Получены турбулентные характеристики как для средних, так и для пульсационных составляющих течения. Результаты получены на основе метода крупных вихрей, они максимально полно описывают реальный физических процесс.
Таким образом, полученные на основе метода крупных вихрей результаты частично могут заменить дорогостоящий реальный физический эксперимент.
Список литературы
fl] Ferziger J.H. Large eddy simulation of turbulent flows // AIAA J. 1977. Vol. 15, N 9. P. 1261-1267.
[2] Sagaut P. Large eddy simulation for incompressible flows. Heidelberg: Springer-Verl., 2002. 423 p.
[3] Яненко H.H. Метод дробных шагов решения многомерных задач математической физики. Новосибирск: Наука, 1967. 196 с.
[4] Самарский А.А., Николаев Е.С. Методы решения сеточных уравнений. М.: Наука, 1978. 592 с.
[5] Zang Y., Street R.L., Koseff J.R. A dynamic mixed subgrid-scale model and its application to turbulent recirculating flows // Phys. Fluids. 1993. Vol. A 5 (12), Dec. P. 3186-3195.
[6] Blay D., Mergui S., Niculae C. Confined turbulent mixed convection in the presence of a horizontal buoyant wall jet // ASME HTD. 1992. Vol. 213: Fundamentals of Mixed Convection. P. 65-72.
Поступила в редакцию 24 апреля 2007 г., в переработанном виде — 21 сентября 2007 г.
