Спектральный анализ вихреразрешающих моделей турбулентности Текст научной статьи по специальности «Физика»

УДК 532.517.4+519.6

СПЕКТРАЛЬНЫЙ АНАЛИЗ ВИХРЕРАЗРЕШАЮЩИХ МОДЕЛЕЙ ТУРБУЛЕНТНОСТИ

ШАКЛЕИН А. А., КАРПОВ А. И., АЛЬЕС М. Ю.

Удмуртский федеральный исследовательский центр Уральского отделения РАН, 426067, г. Ижевск, ул. Т. Барамзиной, 34

АННОТАЦИЯ. Работа посвящена исследованию моделей турбулентности с позиции их возможностей адекватного разрешения широкого спектра частот вихревого течения, что, в свою очередь, является необходимым условием при решении комплексных задач, сочетающих разнообразные физико-химические процессы, например, при горении. Основной акцент в работе делается на модель k-ш SST DDES, хорошо зарекомендовавшую себя в области вихреразрешающего моделирования пристенных турбулентных течений. Верификация осуществляется на базе экспериментальных данных, полученных Comte-Bellot и Corrsin (CBC), на задаче затухания изотропной турбулентности. Представлены спектры турбулентной энергии для различных моделей и аппроксимационных схем. Показано, что модель k-ш SST DDES адекватно разрешает турбулентное течение.

КЛЮЧЕВЫЕ СЛОВА: турбулентное течение, изотропная турбулентность, энергетический спектр, вихреразрешающее моделирование.

ВВЕДЕНИЕ

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

Для повышения результативности проводимых комплексных исследований, связанных с увеличением энергоэффективности установок, предсказанием характера протекания и последствий катастроф, снижения горючести строительных материалов, необходимо адекватно описывать все протекающие в рамках рассматриваемой задачи процессы, в том числе, турбулентное течение [1, 2]. Таким образом, следует применять ориентированные на вихревое течение модели, которые позволяют точно описать основные характерные эффекты: анизотропию крупных вихрей, каскадный механизм переноса энергии, диссипацию мелкомасштабных пульсаций, а также локальные характеристики на стенке (трение и теплообмен) и переходные режимы течения.

Стоит отметить, что модели DES, в частности основанные на модели k-ю SST, ориентированы на подобный класс задач и позволяют как получить информацию о крупных энергосодержащих вихрях, так и без высоких вычислительных затрат в рамках RANS подхода описать пристеночную область течения.

Таким образом, данное исследование направлено на изучение вихреразрешающих моделей с точки зрения оценки способности моделей и используемых методов адекватно разрешать весь возможный спектр вихревого течения [3 - 5]. В работе рассматриваются следующие модели: quasi-DNS (qDNS) [6], k-ю SST DDES [7], модель Смагоринского и модель LES с одним уравнением для кинетической энергии подсеточной турбулентности (k LES) [8]. Проверка моделей осуществляется на основе сравнения результатов численных

расчетов с экспериментом [9] по спектральному распределению характеристик турбулентности. Используемый в работе подход позволяет оценить, насколько хорошо разрешаются моделью различные масштабы турбулентного течения.

МАТЕМАТИЧЕСКАЯ МОДЕЛЬ

В соответствии с методом qDNS течение описывается системой уравнений Навье-Стокса

Эи, Эх,-

= 0.

Эй,- Эй,-—1- + и,-—1-Эt J Эх

:-1 ЭР + .

Э Эи,-

V-

Эх

(1) (2)

V Pdxi dxj -j

Вихреразрешающее моделирование проводится в рамках подходов LES и DES, основанных на фильтрации полных уравнений Навье-Стокса (1) - ( 2)

^ (3)

Эх,-Эи,

Эt

+ и

Эи-

Эх

1 Эр Э

Р Эхi Эх J

V+V

)duL

sgs I Эх-

(4)

В соответствии с моделью Смагоринского подсеточная вязкость рассчитывается в виде \2

V

sgs

= (CsA)2

S,;

(5)

где Cs = 0,167.

В выбранной модели k LES уравнение относительно кинетической энергии подсеточной турбулентности выглядит следующим образом

Э ( ^ Эк

-->- S ■■

sgs У

Эк _ Эк

--+ и --=-

Эt Эх,- Эх,-

„ 772 2 Эик , к V + Vsgs)—+ 2vSgSSjк -Ce-

3/2

sgs' Эх-

3 Эхъ

' J J J " к

при этом подсеточная вихревая вязкость рассчитывается как

A

(6)

Vsgs = СкАк

1/2

(7)

где Ck = 0,094 и Ce = 1,048.

Модель турбулентности k-ю SST DDES замыкает исходную систему уравнений (3) - (4) двумя уравнениями [7]

Эк _ Эк Э / ч Эк

- + и J

Эt J Эх-Эю _ Эю

--+ и J-

Эt JЭх,

Эх] Э

(V + StVsgs )

Эх

+ min

G,c1b кю -b кwF,

DDES-

Эх

(V + SwVsgs )

+ <-pV-^ Ю(1 - Fi ).

J J J

Подсеточная вязкость определяется в следующем виде a^i

V =

sgs

max

a1w,

S,,

Fo

(8) (9)

(10)

Функция БООЕ§ в диссипационном слагаемом уравнения переноса кинетической энергии турбулентности имеет вид

F

DDES

= max

I

к-ю

С

DDES

A

(1 - Fi ) ,1

(11)

Константы модели определяются в виде весовой комбинации ф = Бф + (1 - Б )ф2 и имеют следующие значения: = 0,85; ак2 = 1,0; аЮ1 = 0,5; аЮ2 = 0,856; а1 = 0,55;

а2 = 0,44; Ь = 0,075; Р2 = 0,0828. Остальные постоянные имеют следующие величины:

Ь = 0,09; а = 0,31; с = 10,0; = 0,61 [7].

Задача решается для области в виде куба с периодическими граничными условиями:

Ф| х = х =Ф|х =х , где ф = {щ,к,.

7 7 ,шт л -^шах

Начальные условия имеют вид щ | 0 = и0, = , ю|0 = . МЕТОД РЕШЕНИЯ

Порядок проведения исследования выглядит следующим образом. Для численного решения исходных дифференциальных уравнений используется программный пакет ОрепБОАМ [10] с открытым исходным кодом, основанный на методе конечных объемов. Временная производная в исходных уравнениях (Эф / Эt) аппроксимируется неявными

Эф фи+1 -фи Эф 3фи+1 - 4фи + фи-1 схемами по времени первого — =- и второго порядков — = р

dt Dt dt 2 At

Конвективные слагаемые в исходных уравнениях (uj df / dxj) аппроксимируются схемами

первого порядка и TVD схемами "limitedLinear" и "filteredLinear".

При использовании противопоточной схемы неизвестные значения расчетной величины на грани ячейки контрольного объема переносятся из смежной ячейки, расположенной со стороны набегающего потока газовой среды.

Подход TVD позволяет ограничивать нестабильности схемы высокого порядка с помощью смешения с противопоточной схемой в следующем виде

ff =fLO + V[fHO -fLOL где fLo - значение функции, рассчитанное по противопоточной схеме, фно - значение функции, рассчитанное по центрально-разностной схеме второго порядка, а функция y определена в виде y = max[min[lim,l],ü] .

Для схемы "limitedLinear" lim = 2r / к , где к = 1 и

(df / dx). (x-+i - x) >

2-Г^--1, ifu.+1/2 > 0,

ф+1-ф

r = <

2 И/Эх)f+i(- x) _ i ifuг+1/2 < о.

fi+1 -fi

В соответствии со схемой "filteredLinear

lim =

min (max [Afi+1/2 - 2Afi,0] ,max [Afi+i/2 - 2Afi+1,°]) .

I - к---J.-¡-T-- ,if Аф:+1/2 > 0,

max (I Afi+1/21, 2Afi ,2Af,-+1) Vi+1/2

min(max[2Afi -Afi+1/2,0],max[2Afi+1 -Afi+1/2,°])

1 -к-fül-1 ОЛЛ олл \--,if Afi+1/2 < 0,

max (Afi+1/2 ,2Afi ,2Afi+1)

где Дф/+1/2 =фг+1 -фу, ф = (Эф/Эх). (х+1 -x), Af¿+1 = (Эф/Эх).+1 (х+1 - х), а k = 0,2 и I = 0,05.

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

Начальное поле скорости генерируется на основе экспериментального энергетического спектра E (к) по методике, приведенной в [11].

Порядок построения энергетического спектра для изотропной однородной турбулентности для области в виде куба с периодическими граничными условиями следующий [3].

- Отразить компоненты скорости й*(х-) на спектральное пространство й*(к-). Данную

процедуру можно реализовать, например, с помощью алгоритмов быстрого преобразования Фурье.

- Вычислить кинетическую энергию по спектрам по каждому направлению Е(к- ) = 1 (й* (к- )й*(к- ) ^ , где символ * - комплексное сопряженное.

- Определить функцию энергетического спектра Е (к) интегрированием кинетической энергии Е(к-) по сфере радиуса к, где к представляет собой модуль вектора с

{к х, ку, к1}.

компонентами {кх, ку, к2 ,

ИСХОДНЫЕ ДАННЫЕ

Верификация проводится на задаче затухания изотропной турбулентности на основе экспериментальных данных [9]. В работе представлены спектры вихревого течения для трех моментов времени: = 0 с; ^ = 0,285 с; ^ = 0,665 с.

В качестве расчетной области выбирается куб с ребром Ь = 0,565 м. Строится декартовая равномерная расчетная сетка с количеством узлов вдоль ребра N = 64. Шаг расчетной сетки (Ь / N) определяется из предположения, что ориентировочно 80 % энергии турбулентности должно разрешаться, в то время как оставшиеся 20 % должны оставаться в подсеточной моделируемой области [3]. Таким образом, на данной сетке разрешается около 76 % кинетической энергии турбулентности и рассеивается 11 % энергии, при этом в подсеточной области моделируется 24 % кинетической энергии турбулентности и рассеивается 89 % энергии соответственно.

По известному экспериментальному спектру для момента времени ^ генерируется

начальное поле скорости. Исходный спектр и спектр синтетической турбулентности представлен на рис. 1. На рисунке (1 и далее) вертикальными пунктирными линиями отмечены границы численного решения к0 =11,1 1/м, что соответствует максимальному размеру расчетной области, и ки = 355 1/м, что соответствует пределу Найквиста, характеризующему минимальный размер вихря, который способна разрешить используемая сетка.

Рис. 1. Спектральное распределение турбулентной энергии в начальный момент времени (эксперимент и сгенерированная на его основе синтетическая турбулентность для проведения расчетов)

Значения к0 и w0, используемые в моделях k-ю SST DDES и k LES, несут в себе информацию о той части спектра турбулентности, которая оказалась отброшена вследствие использования расчетной сетки с шагом большим Колмогоровского масштаба длины. Для определения этих величин необходимо проинтегрировать исходный экспериментальный спектр от предела Найквиста ( кп ) до волнового числа, соответствующего Колмогоровскому

масштабу длины ( кл ) к0 = Гк E (к) dk и e = Г к 2vk E(к)dk. В результате получаем

к0 = 0,0188 м2/с2 и e0

2 3 e 0

0,42 м/с . Учитывая, что w = —, получаем w

Ь к

:247 1/с.

РЕЗУЛЬТАТЫ

Проведенные расчеты показывают (рис. 2), что для моделей LES Смагоринского и qDNS происходит заметное перераспределение энергии: в инерциальном диапазоне частот интенсивность вихрей, отвечающих за передачу энергии от крупных к мелким в соответствии с механизмом каскадного переноса энергии, занижается, в то время как в области высоких волновых чисел энергия завышена, что объясняется тем, что модели вносят недостаточное количество численной модельной вязкости. Модели k-ю SST DDES и k LES показывают хорошее согласование с экспериментом.

На рис. 3 приведены спектры, полученные из расчетов с использованием различных схем по модели k-ю SST DDES. Так, результат, полученный с использованием схемы первого порядка аппроксимации по времени (кривая 2, рис. 3) отличается от эксперимента и расчета по схеме второго порядка в области мелкомасштабной турбулентности, что свидетельствует о внесении дополнительной численной вязкости схемой первого порядка. Противопоточная схема и схема TVD "limitedLinear" создают большую искусственную вязкость, что существенным образом искажает результаты решения не только в мелкомасштабном диапазоне, но и в интегральной области спектра.

Рис. 2. Спектральное распределение турбулентной энергии в момент времени t2 для различных моделей: 1 - k-ю SST DDES; 2 - k LES; 3 - LES Смагоринского; 4 - qDNS

10-;

10-'

ю-

10 '

101

......... ......... : ехр - : 1 -------- : " 2 .............. 3......

1/ 1 \ \ -NV. -i \ \ ' \ Ч

- \ \ I \ > -! \ \ \ ........| 1 X \\\ \\ \ \\ \ \\ \

102 к, [1/m]

103

Рис. 3. Спектральное распределение турбулентной энергии в момент времени t2 для различных аппроксимационных схем:

1 - второй порядок по времени, "filteredLinear";

2 - первый порядок по времени, "filteredLinear";

3 - второй порядок по времени, "limitedLinear";

4 - второй порядок по времени, upwind

Дальнейшие расчеты проводятся по модели k-ю SST DDES с использованием схемы второго порядка по времени и TVD схемы "filteredLinear" для конвективных слагаемых. Применение данных схем показывает хорошее согласование получаемых результатов

с экспериментом как по спектральному распределению энергии и скорости диссипации, так и кинетической энергии разрешаемой турбулентности (рис. 4 - 6).

ю2 ю3

к, [1/т]

Рис. 4. Изменение спектральной турбулентной энергии во времени (эксперимент и расчет)

ю2 ю3

к, [1/т]

Рис. 5. Изменение спектральной скорости диссипации турбулентной энергии во времени (эксперимент и расчет)

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 I, [5]

Рис. 6. Изменение турбулентной энергии во времени (эксперимент и расчет)

На рис. 7 представлено распределение модуля скорости в расчетной области в различные моменты времени. По представленным полям скорости наблюдается затухание вихревого течения.

Рис. 7. Поля модуля скорости в различные моменты времени: а - Ь - в -

ЗАКЛЮЧЕНИЕ

Рассматриваемый в работе режим течения является отчасти идеальным в том плане, что течение с изотропной однородной турбулентностью редко встречается в практических задачах из-за преобладания анизотропных эффектов вследствие взаимодействия потоков с твердыми границами. Также следует отметить, что турбулентное течение здесь подчиняется второй гипотезе подобия Колмогорова, в рамках которой можно четко выделить инерциальный диапазон энергетического спектра, что позволяет рассматривать данный подход с точки зрения проверки возможности отражения моделью каскадного механизма передачи энергии. Однако, при различных особенностях течений (например, низкие числа Рейнольдса, переходный режим течения) общая картина будет существенным образом отличаться от универсальной, что, в первую очередь, затруднит проведение спектрального анализа турбулентности. Несмотря на вышесказанное, проведенное исследование позволило оценить способность рассматриваемых моделей (LES с одним уравнением для k и k-ю SST DDES) адекватно разрешать вихри различных размеров, что в конечном итоге позволит применять модель для описания тонких эффектов влияния турбулентного течения на различные физико-химические процессы.

Исследование выполнено при финансовой поддержке РФФИ (проект № 16-08-00110 а). СПИСОК ЛИТЕРАТУРЫ

1. Novozhilov V., Joseph P., Ishiko P., Shimada T., Wang H., Liu J. Polymer combustion as a basis for hybrid propulsion: a comprehensive review and new numerical approaches // Energies, 2011, vol. 4, pp. 1779-1839.

2. Ren N., Wang Y., Vilfayeau S. Trouve A. Large eddy simulation of turbulent vertical wall fires supplied with gaseous fuel through porous burners // Combustion and Flame, 2016, vol. 169, pp. 194-208.

3. Pope S. B. Turbulent flows. Cambridge University Press, 2000. 773 p.

4. Rozema W., Bae H. J., Moin P., Verstappen R. Minimum-dissipation models for large-eddy simulation // Physics of Fluids, 2015, vol. 27, no. 085107.

5. Vilfayeau S. Large eddy simulation of fire extinction phenomen: diss. ... PhD. University of Maryland, 2015.

309 p.

6. Komen E. M. J., Camilo L. H., Shams A., Geurts B. J., Koren B. A quantification method for numerical dissipation in quasi-DNS and under-resolved DNS, and effects of numerical dissipation in quasi-DNS and under-resolved DNS of turbulent channel flows // Journal of Computational Physics, 2017, vol. 345, pp. 565-595.

7. Menter F. R. Kuntz M., Langtry R. Ten years of industrial experience with the SST turbulence model // Turbulence, Heat and Mass Transfer, 2003, vol. 4, pp. 625-632.

8. Menon S., Yeung P.-K., Kim W.-W. Effect of subgrid models on the computed interscale energy transfer in isotropic turbulence // Computers & Fluids, 1996, vol. 25, no. 2, pp. 65-180.

9. Comte-Bellot G., Corrsin S. Simple Eulerial time correlation of full- and narrow-band velocity signals in grid-generated, 'isotropic' turbulence // Journal of Fluid Mechanics, 1971, vol. 48, no. 2, pp. 273-337.

10. Weller H. G., Tabor G., Jasak H., Fureby C. A tensorial approach to computational continuum mechanics using object-oriented techniques // Computers in physics, 1998, vol. 12, no. 6, pp. 620-631.

11. Saad T., Cline D., Stoll R., Sutherland J. C. Scalable Tools for Generating Synthetic Isotropic Turbulence with Arbitrary Spectra // AIAA Journal, 2017, vol. 55, no. 1, pp. 327-331.

SPECTRAL ANALYSIS OF LARGE EDDY TURBULENT MODELS

Shaklein A. A., Karpov A. I., Alies M. Yu.

Udmurt Federal Research Center, Ural Branch of the Russian Academy of Sciences, Izhevsk, Russia

SUMMARY. Many flows observed in practical applications have turbulent character with different eddies lengthscales presented and interacted with each other as well as with enclosures. These interactions result in complex effects -intensification of transfer due to turbulent mixing, large eddy anisotropy - having influence on compound physical and chemical processes. In combustion applications (e.g. fires, power plants, filtration combustion) sufficient amount of processes have to be considered: heat and mass transfer, turbulent combustion, radiative energy transfer, conjugate interaction between gas phase and solid bodies. In order to improve predictions of combustion applications parameters

primary processes have to be resolved properly, especially turbulent flow. Thus, the study is aimed on investigation of turbulent models from the point of view of proper resolving of wide range of turbulent frequencies. In the study several models are examined, namely, qDNS, LES with transport equation for kinetic turbulent energy, LES Smagorinsky, with k-ffl SST DDES being the primary objective. The mathematical formulation of models considered is presented in details. Boundary and initial conditions are given. Numerical algorithm for solving governing differential equations is based on open source OpenFOAM software formulated on the finite volume method. Since numerical schemes are known to have substantial impact on results, especially for large eddy simulations, temporal and convective terms approximation procedure was examined, with diffusion term being fixed to be approximated with linear scheme. Thus, first (Euler) and second (backward) order fully implicit time schemes and first (upwind) and second (bounded, TVD) order convection schemes were applied. The initial velocity field was generated based on existing energy spectrum collected by Comte-Bellot and Corrsin (CBC) from experiment. An algorithm for turbulent energy spectrum calculation from computed velocity field is given in details. Spectral analysis presented shows several flaws of qDNS and LES Smagorinsky models resolving turbulent fluctuations, namely, lowering eddy energy in inertial subrange and overestimating eddy energy of small lengthscales because of insufficient artificial viscosity predicted. On the contrary, LES with one equation for k and k-m SST DDES models predicted turbulent spectrum very accurately in comparison with experimental data. It has to be noted, that considered turbulence is somewhat ideal from the point of view of being isotropic. Indeed, most of engineering application flows have anisotropic structure because of interactions with walls. Nevertheless, carried out study showed ability of models to properly resolve wide range of turbulent frequencies and their applicability for simulation of complex physical and chemical processes interactions.

KEYWORDS: turbulent flow, isotropic turbulence, energy spectrum, large eddy simulation.

REFERENCES

1. Novozhilov V., Joseph P., Ishiko P., Shimada T., Wang H., Liu J. Polymer combustion as a basis for hybrid propulsion: a comprehensive review and new numerical approaches. Energies, 2011, vol. 4, pp. 1779-1839. doi: 10.3390/en4101779

2. Ren N., Wang Y., Vilfayeau S. Trouvé A. Large eddy simulation of turbulent vertical wall fires supplied with gaseous fuel through porous burners. Combustion and Flame, 2016, vol. 169, pp. 194-208. doi: 10.1016/j.combustflame.2015.12.008

3. Pope S. B. Turbulent flows. Cambridge University Press, 2000. 773 p.

4. Rozema W., Bae H. J., Moin P., Verstappen R. Minimum-dissipation models for large-eddy simulation. Physics of Fluids, 2015, vol. 27(8), no. 085107. doi: 10.1063/1.4928700

5. Vilfayeau S. Large eddy simulation of fire extinctionphenomen. University of Maryland, 2015. 309 p.

6. Komen E. M. J., Camilo L. H., Shams A., Geurts B. J., Koren B. A quantification method for numerical dissipation in quasi-DNS and under-resolved DNS, and effects of numerical dissipation in quasi-DNS and under-resolved DNS of turbulent channel flows. Journal of Computational Physics, 2017, vol. 345, pp. 565-595. doi: 10.1016/j.jcp.2017.05.030

7. Menter F. R. Kuntz M., Langtry R. Ten years of industrial experience with the SST turbulence model. Turbulence, Heat and Mass Transfer, 2003, vol. 4, pp. 625-632.

8. Menon S., Yeung P.-K., Kim W.-W. Effect of subgrid models on the computed interscale energy transfer in isotropic turbulence. Computers & Fluids, 1996, vol. 25, no. 2, pp. 65-180. doi: 10.1016/0045-7930(95)00036-4

9. Comte-Bellot G., Corrsin S. Simple Eulerial time correlation of full- and narrow-band velocity signals in grid-generated, 'isotropic' turbulence. Journal of Fluid Mechanics, 1971, vol. 48, no. 2, pp. 273-337. doi: 10.1017/S0022112071001599

10. Weller H. G., Tabor G., Jasak H., Fureby C. A tensorial approach to computational continuum mechanics using object-oriented techniques. Computers in physics, 1998, vol. 12, no. 6, pp. 620-631. doi: 10.1063/1.168744

11. Saad T., Cline D., Stoll R., Sutherland J. C. Scalable Tools for Generating Synthetic Isotropic Turbulence with Arbitrary Spectr. AIAA Journal, 2017, vol. 55, no. 1, pp. 327-331. doi: 10.2514/1.J055230

Шаклеин Артем Андреевич, кандидат физико-математических наук, научный сотрудник Института механики УдмФИЦ УрО РАН, e-mail: mx.oryx@,gmail.com

Карпов Александр Иванович, доктор физико-математических наук, главный научный сотрудник Института механики УдмФИЦ УрО РАН, e-mail: karpov@udman. ru

Альес Михаил Юрьевич, доктор физико-математических наук, профессор, директор УдмФИЦ УрО РАН, e-mail: aliesmy@mail. ru