CC BY
10
УДК 517.958:531.35, 519.6, 662.612.3
ПРЯМОЕ ЧИСЛЕННОЕ МОДЕЛИРОВАНИЕ ТУРБУЛЕНТНЫХ ТЕЧЕНИЙ И ГОРЕНИЯ С ИСПОЛЬЗОВАНИЕМ РАЗНОСТНЫХ СХЕМ ВЫСОКОГО ПОРЯДКА ТОЧНОСТИ
КИСАРОВА С.Ю. , КОРОЛЕВА М.Р., ГАБДУЛЛИН Р.Р.
Институт прикладной механики УрО РАН, Ижевск, Россия, prim'audman. ru
АННОТАЦИЯ. Турбулентные течения являются объектами многих исследований. В настоящее время для расчета турбулентных потоков разрабатываются разностные схемы высокого порядка точности, позволяющие детально описать особенности таких течений. В работе представлен класс разностных схем высокого порядка для прямого численного моделирования турбулентных течений и исследования процессов горения однородной смеси в каналах различной конфигурации.
Для прямого численного моделирования нестационарных турбулентных течений сжимаемого газа использовалась система полных нестационарных, трехмерных уравнений Навье-Стокса. Для исследования течений инертного и реагирующего газа -система уравнений Навье-Стокса для химически реагирующей смеси:
dQ dF dG dH
- + - + - + -:
dt dx dy dz
(
dD dB dA
—+—+ —
dx dy dz
\
+B
(1)
у
Вектора Q, F, G, H, D, B, A, Bi записываются следующим образом:
Q
Г p^
pu pv pw
pE
p у
f
F =
pu
.2
Л
Г о ^
D
xy
V Dx у
pu + p puv puw puH
Р11Сг у
Г о ^
G =
B
yy
Чу
D у
pv puv pv2 + р pvw pvH
v PvC У
Г о ^
H
pw puw
pvw
2
A =
4z
V D у
Bi
pw + p pwH PwCi у
Г о ^ о о о
Bh
V- BC у
Eak/
Bc =P C1C2 Ke , BH = QeffBc.
В системе уравнений (1) x, y, z - декартовы координаты; p - плотность газа; и, v, w - компоненты вектора скорости газа; T - температура газа; p - давление газа; Е, H -полные энергия и энтальпия газа; K, Eak, Qeff - предэкспонент, энергия активации и тепловой эффект газофазных реакций; тхх, Txy, Tz, %, Tyz, Т - компоненты тензора 362
т
т
т
xx
xy
xz
т
т
т
т
xz
zz
Ч
x
Г тЛ05
вязких напряжений; qx, qy, qz - тепловые потоки в соответствующих координатных направлениях: Сг - концентрация /-го компонента в газовой смеси; йх, Ьу, 02 -диффузионные потоки концентрации.
Зависимость величины динамической вязкости от температуры учитывалась следующим образом:
' Т
Т VТ0 у
Система уравнений (1) записана в безразмерных переменных. Все основные геометрические параметры отнесены к входной ширине канала Л, а независимые переменные - к соответствующим параметрам потока на входе.
Интегрирование системы (1) по времени проводилось методом Рунге-Кутта третьего порядка точности:
б0 = б",
61 = б" -Д*ф"),
62 = 4 б" + 4 б1 - 4 Ф1), (2) б3=3 б"+3 б2 - 2 Ф2),
Q = Qn+\
где L(Q) - оператор, включающий пространственные производные системы уравнений (1). Использование данной схемы позволило значительно увеличить временной шаг, что ускорило процесс получения результирующих данных.
Для аппроксимации пространственных производных использовалась комбинация разностных схем высокого порядка точности [1].
Рассмотрим, как проводился расчет пространственных производных с заданным порядком точности N на примере потока F исходной системы (1). Для аппроксимации потоков использовался метод конечных объемов:
f—1 = F+1/2 - F-1/2 (3) ^ dx J Ax '
Вычисление потоков на гранях ячеек осуществлялось следующим образом:
- Q(-)
где
а = max
/ dQ
F'+1 / 2 = f (Qi+1/2 ) ,
N/2
F+1/2 = F+1/2 -1 а(qI+1/2 - Q^/2), (4)
Qi+1 / 2 = Li+1 / 2 (Qi ) = Z (Qi+m + Qi-m+1 ) .
1/2 _ ^¿+1/2\^г/ = Е
т=1
Величины б(+) и б(-) выражения (4) на гранях вычисляются согласно алгоритму
WENO по следующим формулам с порядком точности к :
к-1 к-1
б/У/2Е с^о/-г+;, (6)
г=0 ] = 0
бй/2 = ЕХ ЕЕ сд-г-;. (7)
г=0 ;=0
Порядок точности к в данном случае выбирался таким образом, чтобы не понизить порядок точности N исходной схемы.
Такая разностная схема обладает достаточно небольшой схемной вязкостью, величина которой не оказывает существенного влияния на результаты численных расчетов и автоматически обеспечивает гладкость решения.
Алгоритм тестировался на одномерной задаче о распаде произвольного разрыва, применялся для прямого численного моделирования турбулентных течений в прямолинейных каналах, каналах с резким расширением и криволинейных каналах (диффузорах), а также использовался в расчетах горения однородной смеси в канале с обратной ступенькой [2].
Во всех проводимых численных экспериментах граничные условия задавались следующим образом. На твердых стенках использовались условия прилипания и непроницаемости, на боковых границах - условия периодичности. На входе в канал для основных параметров потока задавались соотношения, позволяющие контролировать расход газа на входе в канал, учитывать влияние возмущений на параметры потока, а также компенсируют эффекты генерации звуковых волн на выходной границе при пересечении ее вихрями [3].
Начальными данными являлось установившееся течение газа (инертного газа при исследовании процессов горения) в плоском канале, полученное в двумерном случае при низком числе Рейнольдса (порядка 200). Поля основных параметров потока распространялись на всю область интегрирования по ширине области интегрирования. Поскольку при двумерных начальных и граничных условиях поток остается двумерным, то для того чтобы течение в двумерной расчетной области стало трехмерным, в поток вносились возмущения. На входной профиль скорости накладывались трехмерные случайные возмущения:
и/ = («,■) + у\и ¿|, /=1> 2 3
где и - вектор скорости, иг - / компонента вектора скорости, (иг) - исходный профиль скорости, |и¿| - модуль вектора скорости, I - интенсивность пульсаций, у -случайна величина, подчиняющаяся нормальному закону распределения, согласно которому у = 0 и л[у = 1.
Возмущения вносились в течение некоторого количества итераций (100 итераций), а число Рейнольдса увеличивалось до нужного значения. Расчет трехмерного течения выполнялся до выхода на стационарный (в осредненном смысле) режим. Затем, при моделировании турбулентных течений с химическими реакциями, на короткое время включался режим подогрева смеси и производился расчет течения с горением. При этом рассчитывались как осредненные, так и пульсационные характеристики потока. Осреднение проводилось в течение 100000 итераций.
С помощью представленного алгоритма расчета турбулентных потоков было проведено детальное исследование течения газа в несимметричном диффузоре. Осредненная картина течения в таком включает в себя следующие характерные зоны: 1) оторвавшийся сдвиговой слой; 2) область рециркуляции; 3) область присоединения сдвигового слоя к стенке; 4) область релаксации.
Было исследовано влияние угла раскрытия диффузора на положение и размер отрывной зоны. При угле раскрытия а = 2° осредненное течение в диффузоре является безотрывным - рис. 1а. При угле раскрытия а = 5° возникает отрывное течение -рис. 1б. Зона отрыва располагается вдоль верхней наклонной стенки диффузора. При дальнейшем увеличении угла раскрытия (а = 10°) длина зоны отрыва существенно не меняется, а сама рециркуляционная область смещается к горловине диффузора -рис. 1в. При этом вихрь в поперечном направлении увеличивается примерно в полтора раза.
Исследовалось также влияние различных граничных условий на картину течения в несимметричном диффузоре с углом раскрытия а = 10°. Расчеты проводились при ширине расчетной области 8И. При этом диффузор в данном направлении ограничивался либо непроницаемыми стенкам, либо условиями периодичности. В первом случае, когда область интегрирования была ограничена непроницаемыми стенками, картина течения не соответствовала экспериментальной. Осредненная картина течения являлась трехмерной и характеризовалась наличием продольных пространственных вихревых образований в углах расчетной области, располагающихся симметрично относительно центра диффузора. При использовании периодических граничных условий осредненная картина течения в диффузоре имеет вид, показанный на рис. 1в. Картина течения в данном случае соответствует экспериментальной [4]. Присутствует зона отрыва, располагающаяся вдоль отклоненной стенки, а также повторное присоединение потока к верхней стенке выходного канала диффузора.
Построенные коэффициенты восстановления давления и сопротивления на наклонной стенке диффузора, осредненные и пульсационные профили скорости также соответствуют эксперименту. На рис. 2 представлены экспериментальные и расчетные
профили средней продольной компоненты скорости в различных точках по оси х. Видно, что в расширяющейся части канала и в начале выходного канала экспериментальные и расчетные данные хорошо согласуются между собой. Незначительное расхождение данных наблюдается лишь вблизи выходного сечения диффузора.
На рис. 3 - 4 приведены графики пульсаций продольной компоненты скорости потока и напряжений Рейнольдса в различных сечениях в первой половине диффузора. Профили среднеквадратичных пульсаций продольной компоненты скорости имеют два ярко выраженных максимума, положение которых определяется наличием возвратного течения вдоль верхней наклонной стенки диффузора и присутствием вихрей в отрывающемся сдвиговом слое.
Рис. 2. Профили средней скорости X +10 • и / V: о - эксперимент,-расчет
4 -
Рис. 3. Профили пульсаций скорости x +150 • и'u'/V: о - эксперимент,-расчет
Рис. 4. Профили напряжений Рейнольдса x + 750 • u' у' / V: о - эксперимент,-расчет
Профили касательных напряжений Рейнольдса в разных сечениях после расширения аналогично профилям пульсаций продольной компоненты скорости имеют два ярко выраженных экстремума, положение которых также соответствует расположениям сдвигового слоя и рециркуляционной области. Видно хорошее соответствие расчетных данных экспериментальным.
Разработанная схема применялась также для моделирования течения инертного и реагирующего газа в канале за уступом. При исследовании течения реагирующего газа предполагалось, что протекает одноступенчатая обобщенная реакция, скорость которой вычислялась по уравнению Аррениуса. В расчетах были приняты следующие величины констант реакции: К = 1,94-106, Еак = 20, Qef = 3; 5; 6.
На рис. 5 приведен график изменения давления по времени для течения инертного газа и во время поджига. Хорошо виден период установления течения и его выход на стационарный режим. Резкое изменение давления в момент времени I = 250 соответствует включению режима подогрева смеси, после чего происходит горение газа с величиной теплового эффекта 3. Как видно, характер изменения давления становится нерегулярным. Амплитуда и частота колебаний давления при горении газа возрастают в 1,5-2 раза.
Рис. 5. Изменение давления по времени для течения инертного газа и во время поджига
На рис. 6. - 7. показаны кривые изменения давления по времени при горении с различными тепловыми эффектами. Частота и амплитуда колебаний меняются при различных режимах. Амплитуда колебаний давления при горении газа с тепловым эффектом 5 почти в 3 раза превосходит амплитуду колебаний горения с тепловым эффектом 3. Частота колебаний также значительно возрастает.
Рис. 6. Изменение давления во время горения Рис. 7. Изменение давления во время горения с с тепловым эффектом 3 и 5 тепловым эффектом 5
Более подробное исследование процесса горения в канале показало, что горение в основном происходит внутри вихрей, что соответствует экспериментальным результатам [5], кроме этого процесс горения подавляет процесс вихреобразования.
При дальнейшем увеличении теплового эффекта до 6 характер течения меняется еще сильнее. Образование вихрей в этом случае происходит в основном в головной части канала. Время от времени эта зона размыкает, оказывая сильное влияние на остальной поток газа. Происходит периодическое выгорание всей смеси в канале, которое сопровождается резким повышением давления и скорости истечения газа -рис.8.
Рис.8. Изменение по времени картины течения горящего газа с тепловым эффектом 6
Представленный алгоритм расчета хорошо описывает такие явления, как ударные волны, контактные разрывы, волны разряжения и отраженные волны . Позволяет проводить моделирование турбулентных течений в областях сложной формы, а также учитывать протекающие в них химические реакции и может быть использован для исследования широкого круга сложных течений.
Работа выполнена при финансовой поддержке РФФИ. Проект № 07-01-96074.
СПИСОК ЛИТЕРАТУРЫ
1. Кисаров Ю.Ф., Кисарова С.Ю., Королева М.Р., Шмыкова В.Н. Численное моделирование дозвуковых и сверхзвуковых течений с использованием разностных схем высокого порядка точности. // Сб. докладов международной научной конференции "Фундаментальные и прикладные вопросы механики". - Хабаровск, 2003. - Т. 1. - С. 143-157.
2. Кисаров Ю.Ф., Кисарова С.Ю., Королева М.Р. Некоторые особенности турбулентного течения и горения при обтекании обратной ступеньки // Сб. докладов международной конференции "Нелинейные задачи теории устойчивости и турбулентность". - М., 2004.
3. Липанов А.М. и др. Численный эксперимент в классической гидромеханике турбулентных потоков / А.М. Липанов, Ю.Ф. Кисаров, И.Г. Ключников. - Екатеринбург: УрО РАН, 2001. - 161 с.
4. Obi S., Aoki K., Masuda S. Experimental and computational study of turbulent separating flow in an asymmetric plane diffuser // IX Symposium on Turbulent Shear Flows. - Kyoto, Japan. - 1993. - Paper P305-1.
5. Р.В. Питц, Дж.У. Дейли. Горение в турбулентном слое смешения за уступом // Аэрокосмическая техника. - 1984. - Т.2, №7, С.74-82.
SUMMARY. Turbulent streams are objects of many researches. Now for numerical modelings of turbulent flows are often developed high-order methods, allowing in details to describe the features of such streams. In work a class of numerical high-order schemes for direct numerical simulation of turbulent flows and flows with chemical reactions in channels of a various configuration is presented.
