УДК 519.6
С.М. Дмитриев1, О.Л. Крутякова2, А.С. Козелков 1-3, А.А. Куркин1, В.В. Курулин2, Д.А. Уткин2
ПРИМЕНЕНИЕ ПОЛУЭМПИРИЧЕСКИХ МОДЕЛЕЙ ТУРБУЛЕНТНОСТИ ДЛЯ МОДЕЛИРОВАНИЯ ТУРБУЛЕНТНОЙ КОНВЕКЦИИ
Нижегородский государственный технический университет им. Р.Е. Алексеева1 Всероссийский научно-исследовательский институт экспериментальной физики -Российский федеральный ядерный центр2 Национальный исследовательский ядерный университет МИФИ3
Рассматривается возможность использования различных подходов к моделированию турбулентности в условиях развитой конвекции при числах Рэлея высоких порядков. Для исследования выбран ряд промыш-ленно ориентированных задач, имеющих экспериментальные данные. Показано, что при числах Рэлея от 109 до 1017 применение вихреразрешающей LES модели позволяет увеличить точность моделирования естественной конвекции в сравнении с RANS моделями турбулентности. Наиболее ярко данное преимущество проявляется для случаев конвекции Рэлея-Бенера, при вертикальном перепаде температур с образованием обширной конвективной зоны высокой интенсивности. Использование модели рейнольдсовых напряжений EARSM показано для случаев естественно-конвективного течения в областях с наличием двугранных углов и преобладанием вторичных токов. При моделировании менее интенсивного конвективного течения различия в применяемых подходах к моделированию турбулентности менее значительны. Показано, что с увеличением значений чисел Рэлея возрастают погрешности в определении теплогидравлических характеристик. Для более точного их определения целесообразно использовать вихреразрешающие подходы к моделированию турбулентности.
Ключевые слова: высокоразвитое конвективное течение, моделирование турбулентности, высокие числа Рэлея, конвекция Рэлея-Бенера, теплоотдача.
Введение
Течения жидкостей с развитой турбулентной конвекцией, характеризующейся сильным влиянием гравитационных сил, представляют огромный интерес, но до настоящего времени слабо изучены. Необходимость исследования данных течений возникает в таких отраслях научных и технических знаний, как астрофизика, геофизика, геодинамика, атомная энергетика. К важнейшим задачам, требующим учета турбулентной конвекции, относятся крупномасштабные течения в атмосфере и жидких ядрах планет, течения жидкометаллических теплоносителей в реакторных установках, конвекция в условиях микрогравитации в космических аппаратах и топливных баках, а также другие разнообразные течения в технических изделиях. Особый интерес в исследовании турбулентной конвекции проявляется в атомно-энергетической отрасли. Первостепенное значение, наряду с повышением технико-экономических показателей атомных энергетических станций (АЭС), имеет обеспечение безопасности реакторных установок при прохождении тяжелых аварий, а также проведение мероприятий с целью минимизации их последствий [1]. При тяжелых авариях, сопровождающихся потерей теплоносителя, происходит разрушение активной зоны. Вследствие этого расплавленные элементы конструкции реактора перемещаются на днище корпуса реактора, что приводит к образованию высокотемпературного тепловыделяющего расплава в его нижней части, где он должен быть локализован. В число мер безопасности реакторной установки входит обеспечение необходимого отвода тепла от границ устройства локализации расплава, где теплообмен происходит в режиме естественной конвекции.
К экспериментам в атомной энергетике, где исследовался механизм конвективного теплообмена жидкости в полостях с объемным нагревом при числах Рэлея высоких порядков 1014 -1017 относятся СОРО [2-5] и BALI [6, 7]. Сложность экспериментов данного класса те-
© Дмитриев С.М., Крутякова О.Л., Козелков А.С., Куркин А.А., Курулин В.В., Уткин Д.А.
чений актуализирует развитие численного моделирования. Поскольку в ходе исследований подтвержден ярко выраженный турбулентный характер течения, одной из основных проблем численного моделирования является турбулентность. Интерес вызывает вопрос влияния турбулентности на характер течения и величину граничной теплоотдачи. Существующий уровень развития вычислительных мощностей позволяет моделировать конвекцию в приемлемые сроки с помощью осредненных по Рейнольдсу уравнений Навье-Стокса. При этом турбулентность моделируется в рамках RANS подхода [8, 9]. Это накладывает ограничения при описании турбулизации потоков как для естественной, так и для вынужденной конвекции, вследствие наличия гравитационных сил, моделируемых в RANS подходе эмпирическими зависимостями, что, в свою очередь, приводит к существенной погрешности описания процесса в целом. Данная проблема может быть решена использованием прямого численного моделирования (DNS). Однако проводимые в настоящее время численные эксперименты с помощью DNS для задач естественной конвекции ограничиваются достаточно скромными числами Рэлея: максимальное значение 107 [10]. Особый же интерес в изучении турбулентной конвекции представляют высокие числа Рэлея, для технических приложений - это числа порядка 1014-1017 и выше.
В работе [11] отмечается, что стандартные низкорейнольдсовые RANS модели не всегда применимы для решения задач турбулентной конвекции. Численное моделирование эксперимента BALI с целью верификации и калибровки различных моделей турбулентности проводилось в [12]. В работе [13] показано, что применение низкорейнольдсовой AKN модели турбулентности (описание AKN модели приведено в [13-15]) позволяет снизить погрешность распределения теплового потока. Моделирование с использованием стандартной к-е модели турбулентности показывает, что ее применение ведет к заниженному распределению теплового потока.
Целью данной работы является численное моделирование процессов высокоразвитой конвекции тепловыделяющей жидкости при высоких числах Рэлея и исследование применения подходов к описанию турбулентности для свободно конвективных течений.
В первой части работы представлены результаты моделирования свободно конвективного течения между двумя стенками разной температуры в замкнутой кубической полости, которая представляет наиболее распространенный класс задач естественной конвекции. В основной части работы представлены результаты численного исследования задач естественной конвекции при высоких числах Рэлея на примере моделирования экспериментов COPO и BALI. Показаны результаты применения различных подходов к моделированию турбулентности в сравнении с результатами экспериментальных исследований.
Подходы к моделированию турбулентности
Моделирование турбулентности опирается на несколько основных подходов. Наибольшей популярностью в настоящее время пользуются RANS модели турбулентности, однако они не являются универсальными и подходящими для решения широкого круга прикладных задач, поскольку описывают только осредненные характеристики, что налагает определенные требования к их применимости на практике [8, 16]. Наивысший рейтинг применимости из RANS моделей имеет модель SST (Shear Stress Transport) [17]. Однако, как и в большинстве RANS моделей, в SST используется гипотеза Буссинеска о турбулентной вязкости [18], которая справедлива лишь в случае изотропной турбулентности. Для турбулентной конвекции существенное влияние оказывают анизотропные свойства течения: это связано со сложным характером потока и наличием вторичных токов. Для правильного предсказания структуры подобных течений необходимо использовать RSM модели (Reynolds Stress Modelling) или альтернативные EARSM (Explicit Algebraic Reynolds Stress Modelling), которые учитывают влияние на основное течение всех компонент тензора напряжений Рейнольд-са с помощью нелинейных соотношений, и тем самым способны увеличить точность расчета
турбулентных, анизотропных течений [18]. В данной работе в качестве альтернативы модели SST, рассматривается явная алгебраическая EARSM модель, которая по своим качествам не уступает дифференциальной RSM и является менее затратной с вычислительной точки зрения [19].
Метод крупных вихрей LES (Large Eddy Simulation) позволяет получить хорошие результаты как для присоединенных, так и для отрывных течений, существенно уточняя прогнозирование ряда основных физических процессов [20]. Данный подход способен предоставлять детальную информацию о нестационарных полях флуктуаций скорости, температуры и давления, которые, в свою очередь, могут оказывать влияние на величину и характер граничной теплоотдачи. Однако метод LES налагает большие требования на качество дискретных моделей и существенно повышает, по сравнению с RANS, объем необходимых вычислительных ресурсов. Таким образом, для конвективных течений встает вопрос исследования строения турбулентного потока и применение вихреразрешающих моделей турбулентности для установления связи мелкомасштабных колебаний с теплоотдачей в рассматриваемой зоне.
Численные эксперименты
Все расчеты в настоящей работе осуществлены с помощью пакета ЛОГОС - российского программного продукта инженерного анализа, предназначенного для решения сопряженных трехмерных задач конвективного тепломассопереноса, аэродинамики и гидродинамики на параллельных ЭВМ. Пакет программ ЛОГОС успешно прошел верификацию и показал достаточно хорошие результаты на серии различных гидродинамических задач [15, 22, 23], включая расчеты турбулентных и нестационарных течений [19, 20] и геофизические явления [24-26].
Задача о развитой турбулентной конвекции
Моделирование течения жидкости в кубической полости между двух разнотемпера-турных стенок представляет фундаментальный интерес, поскольку течения с наложенным вертикальным перепадом температуры характеризуются наличием разномасштабных вихревых структур и учитывают основные особенности естественно конвективного течения. В работе приведены результаты численного моделирования циркуляции дистиллированной воды в области, подогреваемой снизу и охлаждаемой сверху. Основные исследования данного вопроса выполнены для конвекции в цилиндрических и кубических полостях. Экспериментальные исследования конвекции в кубической полости с вертикальным перепадом температуры проводились в широком диапазоне чисел Рэлея 103-109.
Экспериментальная установка представляет собой кубическую область со стороной Б = 250 мм [10]. Горизонтальные стенки изготовлены из меди и выступают в качестве теплообменников, а вертикальные стенки изготовлены из плексигласа, обеспечивают изотермические граничные условия. С помощью термостатов через теплообменники пропускается тер-мостатирующая жидкость, которая нагревает (охлаждает) теплообменники относительно средней температуры жидкости. В данной работе рассматривается течение, характеризующиеся перепадом температуры в 20 Х, что соответствует числу Рэлея 6,1 х109, которое определяется по формуле (1):
Ка = *' (1)
Х-у
где g - ускорение свободного падения, в - коэффициент температурного расширения, Н - характерный размер, X - коэффициент теплопроводности, V - кинематическая вязкость, ДГ - перепад температур.
В качестве теплоносителя используется вода, свойства которой заданы зависящими от температуры [27]. Моделирование проведено в рамках вихреразрешающего подхода LES и RANS моделей SST и EARSM.
Использовалась изотропная блочно-структурированная сетка, состоящая из 3,375 млн ячеек правильной гексаэдральной формы с размером 0,0016 м. Сеточная модель соответствует требованию применимости стандартных пристеночных функций с параметром y + не выше 5, а также позволяет разрешить мелкомасштабную турбулентность в рамках LES-расчета.
В качестве схемы дискретизации по пространству и в уравнении сохранения импульса для LES-модели использовалась схема BCD с откалиброванной для данной схемы константой Смагоринского. Описание схем представлено в [9, 21]. Шаг по времени подбирался таким образом, чтобы обеспечить значение числа Куранта порядка единицы. Нестационарный расчет проводился до момента времени 4200 с, что соответствует оценочному времени основных характеристик физических полей в эксперименте. Осреднение основных физических полей проводилось с момента времени 100 с.
На рис. 1 показаны расчетные и экспериментальные поля осредненной скорости и распределения спектральной плотности энергии пульсаций крупномасштабной циркуляции
/г 2 , i2\
- Z7 Vx+ Vz ,
в центральном сечении экспериментальной установки Еpuls avg —--- (рис. 2).
2
а) б)
Рис. 1. Распределение осредненной по времени амплитуды скорости:
а - эксперимент; б - расчет
а) б)
Рис. 2. Распределение осредненной по времени энергии турбулентных пульсаций:
а - эксперимент; б - расчет
Максимальное значение скорости на картине экспериментального мгновенного поля достигается вблизи стенок, в противоположных углах области возникают вихри, численное поле распределения спектральной плотности энергии пульсаций представляет похожую картину, описывая локализованные в углах вихри небольшого размера. В эксперименте и в расчете структура течения очень близка. Хорошо видно, что в области возникает крупномасштабный вихрь, соизмеримый с ее размерами. Помимо основного вихря, на средних полях также существуют вторичные вихри меньшего масштаба, локализованные в противоположных углах кубической области.
Для качественной и количественной оценки на рис. 3 представлены профили горизонтальной и вертикальной компонент скорости осредненной по времени. Как видно из рис. 3, профили скорости демонстрируют хорошее совпадение, повторяя не только структуру потока, но и воспроизводя значение скорости в пограничных слоях.
300 ... ^ — — 250 *** к ^ • • Эксперимент -LES
--SST --EARSM
150^
\Ч v 4 л«
0
Ux, мм/сек
а)
10
-10
• • Эксперимент
-LES
--SST // I
--EARSM dt i if/ 1
^jb* - Toó" 150 200 250
/fi'
1 '£'
1 'p
1 '*
X, мм
Рис. 3. Профили средней по времени профиль горизонтальной (а) и вертикальной (б) компоненты скорости при фиксированной координате х = 125 мм и z = 125 мм соответственно
Профили скорости, полученные с помощью модели линейной вихревой вязкости SST, отклоняются от эксперимента. Среднеквадратичное отклонение по данной модели составило порядка 6 %. Различия в расчетных данных могут быть связаны с геометрической особенностью области, характеризующейся большим количеством двугранных углов. Так, в [19] на задаче течения жидкости в трубе квадратного сечения показано, что модель SST, дает большую погрешность, и для подобных задач показано применение модели EARSM. Использование данной модели позволило значительно сократить вычислительную погрешность, среднеквадратичное отклонение результатов EARSM составило около 2,5 %. С помощью LES модели получено наилучшее согласие с экспериментальными данными. Расчетные профили осредненной скорости по LES с высокой точностью воспроизводят общий характер конвективного теплообмена в кубической полости, подогреваемой снизу и охлаждаемой сверху, относительная погрешность результатов составила ~ 1,5 %. Результаты по всем моделям получены на одной и той же расчетной сетке.
Численное решение данной задачи показывает, что моделирование конвективного движения для чисел Рэлея ~109 по всем моделям турбулентности дает приемлемый результат. Однако практически важные течения, вызывающие большой интерес при обосновании безопасности АЭС характеризуются более высокими числами Рэлея, вплоть до 1017. Применимость представленных подходов к моделированию турбулентности для такого класса течений можно оценить на экспериментах COPO и BALI, которые представлены ниже.
Моделирование эксперимента COPO
В эксперименте COPO рассматривается задача развитой турбулентной конвекции водного раствора соли в замкнутой области. Эксперимент COPO имитирует аварийную ситуацию в водо-водяном энергетическом реакторе с эллиптическим днищем (ВВЭР), при которой происходит разрушение активной зоны, что приводит к образованию высокотемпературного тепловыделяющего расплава в нижней части корпуса реактора [1-5]. Для предотвращения дальнейшего разрушения реакторной установки осуществляется теплоотвод, путем охлаждения нижней части корпуса реактора водой, верхняя часть реактора охлаждается за счет радиационного теплообмена, течение внутри области происходит в режиме естественной конвенции.
Общий вид экспериментальной установки показан на рис. 4: она представляет собой тонкий вертикальный срез соответствующий центру нижней части корпуса реактора ВВЭР-440 с эллиптическим днищем. Эксперимент проводился в масштабе 1:2 в slice-геометрии, длина установки составляла 1,77 м, толщина слоя - 0,094 м, большая горизонтальная ось полуэллипса - 0,89 м, вертикальная ось b = ^ a.
Рис. 4. Схема установки COPO
Жидкостью, имитирующей расплав, являлся раствор соли в воде, объемный нагрев которого обеспечивали распределенные по всей области электроды. Вертикальные границы были теплоизолированными. Эллиптическая граница охлаждалась жидким азотом через прослойку теплоизолятора, что приводило к образованию относительно тонкой корки льда и обеспечивало изотермические граничные условия. В одном из режимов верхняя граница экспериментальной установки была теплоизолированной, что аналогично случаю, когда в устройстве локализации расплава вышерасположенный оксидный бассейн покрыт коркой, и тепловыделение не велико. Для остальной серии экспериментов верхняя граница охлаждалась аналогично эллиптической.
Основными характеристиками течения с внутренними источниками тепла являются параметры плотности теплового потока на границе и максимальное превышение температуры жидкости в объеме над температурой границы. Безразмерные числа, характеризующие данные величины - критерий Нуссельта и Рэлея соответственно. Для тепловыделяющего расплава используют модифицированное число Рэлея (2):
= ^^ (2)
к-аТ -V
где g - ускорение свободного падения, в - коэффициент теплового расширения, qv - внут-риобъемное энерговыделение, R - высота полуэллипса, X - коэффициент теплопроводности среды, ат - коэффициент температуропроводности среды, v - кинематическая вязкость. Значение модифицированного числа Рэлея в данном эксперименте составляло ~ 1015.
Показателем эффективной теплоотдачи через границу объема служит число Нуссельта (3):
Nu = . (3)
IAT
В расчете все теплофизические величины заданы зависящими от температуры [27]. Коэффициент температурного расширения был задан зависящим от температуры и определяется соотношением (4):
Р = (0,007114 • T -1,87241445)/(754,5045 +1,8724 • T - 0,00355711-T2). (4)
Моделирование проведено в рамках вихреразрешающего подхода LES и RANS моделей SST и EARSM. Расчетная область соответствовала полной модели экспериментальной установки. Использовалась блочно-структурированная сетка, состоящая из 2,5 млн ячеек гексаэдральной формы, размер ячеек в основной области 0,004 м, со сгущением ко всем стенкам (размер пристеночной ячейки 0,0002 м). Сеточная модель соответствует требованию применимости стандартных пристеночных функций с параметром y+ не выше 5, а также позволяет разрешить мелкомасштабную турбулентность в рамках LES расчета.
В качестве схемы дискретизации по пространству и в уравнении сохранения импульса для LES-модели использовалась схема BCD с откалиброванной для данной схемы константой Смагоринского. Шаг по времени подбирался таким образом, чтобы обеспечить значение числа Куранта порядка единицы. Нестационарный расчет проводился до момента времени 5000 сек., осреднение основных физических полей проводилось с момента времени 1000 сек.
Численное моделирование выполнено для серии экспериментов COPO, в которых присутствовал режим с теплоизолированной верхней стенкой (режим № 1), а также три режима с охлаждаемой верхней стенкой, различным значением энерговыделения и различной высотой области (аспектным соотношением).
Таблица1
Параметры экспериментальных режимов
Режим № Высота области (м) Режим для верхней стенки Объемное энерговыделение (кВт/м2) Значение модифицированного числа Рэлея
1 0,847 адиабатическая 46,5 1,7 e+15
2 0,639 T=273K 153,2 1,36 e+15
3 0,73 T=273K 161,5 2,9 e+15
4 0,847 T=273K 104,3 3,8 e+15
Экспериментально подтверждается, что для режима № 1 с теплоизолированной верхней стенкой характерно стратифицированное распределение. На рис. 5, а представлено мгновенное распределение температуры для режима №1 (модель LES) видно, что жидкость с более высокой температурой (минимальной плотностью) находится вверху области, возмущения в потоке практически отсутствуют, наблюдается четкая температурная стратификация. На рис. 5, б представлено мгновенное распределение безразмерной температуры (модель LES) для экспериментального режима № 4 с охлаждаемой верхней границей. Для данного режима стратифицированное распределение температуры наблюдается только в нижней части расчетной области, большую часть верхней области занимает зона конвекции Рэлея-Бенера, где наблюдаются крупные вихревые структуры.
Для экспериментального режима № 1 на рис. 6, а приведено распределение теплового потока на охлаждаемой эллиптической границе, в сравнении с экспериментальными данными. Эксперимент содержит данные по распределению плотности теплового потока, как по левой части корпуса, так и по правой. Разница данных графиков, по всей видимости, обусловлена чистотой проведения эксперимента. Модели SST и EARSM представляют похожий результат распределения плотности теплового потока, среднеквадратичная погрешность составляет около 6 %, применение EARSM дает небольшое преимущество, уменьшая погрешность на 0,5 %. Среднеквадратичная погрешность при LES моделировании составила около 3,5 %.
а) б)
Рис. 5. Мгновенное поле безразмерной температуры, модель LES:
а - режим № 1; б - режим № 4
а) б)
Рис. 6. Распределение теплового потока на эллиптической границе: а - режима № 1; б - режим № 4
Для экспериментального режима с охлаждаемой верхней границей (рис.6, б) результаты имеют более высокую погрешность. Все модели завышают тепловой поток в верхней части расчетной области. Результаты, полученные с использованием RANS-моделей, соотносятся между собой, давая приблизительно одинаковый результат, среднеквадратичное отклонение для EARSM и SST-моделей составляет 11 % и 13 % соответственно. Кривая плотности теплового потока, полученная с помощью LES-моделирования, лежит ближе к эксперименту, погрешность составила порядка 6 %.
Основным интегральным показателем рассмотренного класса течений является распределение граничной теплоотдачи, характеризуемое числом Нуссельта. По результатам проведенных за последнее время экспериментальных исследований была получена зависи-
мость числа № от модифицированного числа Рэлея (2), зависимость справедлива для Яя ~ 1013-1017 [3]:
Мы = 0,385 • Яа^'233 (6)
На рис. 7 представлены результаты численного моделирования трех экспериментальных режимов, в соответствии с табл. 1, для которых приведена зависимость числа Нуссельта (3) от модифицированного числа Рэлея (2).
2000
Z
1000
500
1.00Е+15 2.00Е+15 4.00Е+15
Rai
Рис. 7. Зависимость числа Нуссельта от модифицированного числа Рэлея (Rai)
Анализируя полученные результаты, можно предположить, что достаточно большое отклонение связано с занижением теплового потока на верхней границе, в то время как на эллиптической границе тепловой поток завышен. Теплосъем с верхней границы менее интенсивен, в сравнении с экспериментальными данными. На эллиптической же границе он завышен в среднем на 30 %. Максимальное отклонение от экспериментальной зависимости достигает 33 %, для LES максимальное отклонение составило 25 %.
В рамках эксперимента COPO рассматривается конвективное течение с числом Рэлея до 1015 [2-5]. Поскольку течения, сопровождающие аварии в реакторных установках характеризуются числами Рэлея еще более высоких порядков, значительный интерес представляет собой изучение конвекции при числах Рэлея ~ 1017. Задачей, перекрывающим область наиболее важную для проблем безопасности ядерных реакторов, является эксперимент BALI.
Моделирование эксперимента BALI
С целью исследования конвекции тепловыделяющей жидкости при тяжелой аварии реактора PWR во второй половине 1990-х гг. была проведена серия экспериментов BALI [6, 7]. В отличие от экспериментов COPO, эксперименты BALI моделируют процессы, происходящие при тяжелой аварии в реакторе со сферическим днищем и при более высоких числах Рэлея. Выполнено четыре серии экспериментов в разных условиях и конфигурациях в масштабе длин 1:1. Общий вид экспериментальной установки показан на рис. 8. Жидкостью, имитирующей расплав, являлся раствор соли в воде, объемный нагрев которого обеспечивали распределенные по всей области электроды. Боковые стенки были теплоизолированными, сферическая граница охлаждалась жидким азотом через прослойку теплоизолятора, что приводило к образованию относительно тонкой корки льда и обеспечивало изотермические граничные условия. В первой серии экспериментов верхняя граница являлась теплоизолированной, во второй серии экспериментов охлаждалась аналогично сферической.
Моделирование проведено в рамках вихреразрешающего подхода LES и RANS-модели EARSM.
-Аналитика
-И- EARSM —SST
режим №2
Расчетная область соответствовала полной модели экспериментальной установки. Сеточная модель содержала 4,2 млн ячеек, размер ячеек в основной области составлял 0,006 м, размер пристеночной ячейки - 0,0001 м и соответстветствовала требованию применимости стандартных пристеночных функций с параметром y + не выше 5.
Настройки расчета в части схемы, шага по времени и физических параметров жидкости аналогичны настройкам при проведении вычислительного эксперимента COPO. Расчет проведен до 6000 с физического времени, осреднение основных физических полей проводилось с момента времени 1000 с.
Численное моделирование выполнено для серии экспериментов BALI, в которых присутствовал режим с теплоизолированной верхней стенкой (режим № 1), а также три режима с охлаждаемой верхней стенкой, различным значением энерговыделения. В табл. 2 представлены параметры режимов.
Л
хладагент
Рис. 8. Схема установки BALI
Таблица 2
Параметры экспериментальных режимов
Режим № Режим для верхней стенки Объемное энерговыделение (кВт/м2) Значение модифицированного числа Рэлея
1 адиабатическая 21,5 3,6 e+16
2 T=273K 5 6,82 e+15
3 T=273K 21,5 3,6 e+16
4 T=273K 42,7 1,07 e+17
Выходные параметры по тепловому потоку и температуре приводятся в безразмерном виде. Из графика зависимости безразмерной температуры по высоте канала (рис. 9) видно, что физическая картина для первого расчетного случая, с теплоизолированной верхней стенкой аналогична задаче COPO. Результаты первой серии численных экспериментов (режим № 1) показывают, что модель турбулентности EARSM с достаточной точностью воспроизводит общий характер конвективного теплообмена в экспериментальной установке BALI. Распределение безразмерной температуры на вертикальной стенке, полученной в ходе расчета, согласуется с экспериментальными данными для LES-модели турбулентности в большей степени, нежели чем для EARSM (рис. 9, а). Распределение безразмерной плотности тепло-
вого потока на охлаждаемой сферической границе хорошо согласуется с опытными данными в верхней части экспериментальной установки, в ее нижней части наблюдается занижение его величины для всех моделей турбулентности (рис. 9, б).
Среднеквадратичная погрешность при измерении теплового потока составляет около 5 % (EARSM), применяется вихреразрешающая модель, кривая плотности теплового потока при LES моделировании лучше соотносится с экспериментальными данными, среднеквадратичная погрешность результатов составляет около 4 %.
а) б)
Рис. 9. Распределение осредненной по времени безразмерной температуры на вертикальной границе (а) и осредненной по времени безразмерной плотности теплового потока на сферической границе (б), режим № 1
Ниже представлены результаты численного моделирования задачи BALI для режима с охлаждаемой верхней границей (режим № 3). При экспериментальном исследовании наблюдается несколько зон течения: верхний неустойчивый слой однородной температуры, где большую часть занимает зона конвекции Рэлея-Бенера; нижняя зона со стратификацией по температуре, при которой жидкость поднимается в центр с низкой скоростью; пограничный слой, где жидкость в контакте со льдом охлаждается и возвращается с высокой скоростью на дно бассейна. Аналогичная картина наблюдается и при численном моделировании, большую часть области занимает зона крупномасштабной конвекции, где локализуются крупные вихревые структуры, теплообмен происходит в режиме высокоразвитой турбулентной конвекции. Для данного режима на рис. 10 представлены мгновенное распределение безразмерной температуры, режим № 3 для EARSM и LES.
T/Tmax T/Tmax
а) б)
Рис. 10. Мгновенное поле температуры: модель LES (а) и EARSM (б)
Стратифицированное распределение температуры наблюдается только в нижней части расчетной области, большую часть верхней области занимает зона конвекции Рэлея-Бенера и она значительно больше, чем зона, наблюдавшаяся при моделировании эксперимента COPO: это может быть связано как с геометрией установки, так и с увеличенным значением числа Рэлея.
Эксперимент содержит данные по распределению температуры и теплового потока для режима течения с двумя охлаждающимися стенками. На рис. 11 приведены результаты по осредненной плотности теплового потока и распределению безразмерной температуры. Анализ результатов расчета показывает, что модель EARSM занижает величину теплового потока в верхней части расчетной области, результаты LES лежат ближе к эксперименту. Среднеквадратичное отклонение результатов EARSM составляет 11 %, LES - 6%.
На рис. 12 представлены результаты численного моделирования трех экспериментальных режимов с охлаждаемой верхней стенкой для режимов № 2, № 3, № 4 в сравнении с аналитической зависимостью числа Нуссельта от модифицированного числа Рэлея (модели турбулентности EARSM и LES).
а) б)
Рис. 11. Распределение осредненной по времени безразмерной температуры по вертикальной границе (а) и распределение осредненной по времени безразмерной плотности теплового потока на сферической границе (б)
Рис. 12. Зависимости числа Нуссельта от модифицированного числа Рэлея
Максимальная погрешность для числа Нуссельта на верхней границе составила около 35 % (EARSM), результаты с использованием вихреразрешающего подхода лежат ближе
к аналитическому решению, в случае максимальных значений числа Рэлея погрешность составляет 33 %. Численное моделирование качественно соответствует физической картине естественной конвекции. Вследствие высокой разницы температур у верхней стенки возникает крупные вихревые потоки. В нижней части бассейна проявляется устойчивая зона температурной стратификации. Данная картина соответствует течению в экспериментальной установке. Распределение безразмерной температуры на вертикальной границе описывается с достаточной точностью, при этом наибольшие отклонения наблюдаются в центральной области. Расчетное распределение плотности безразмерного теплового потока на охлаждаемой верхней границе описывается с более значительным отклонением от экспериментальных данных.
Заключение
Численное моделирование процессов турбулентной естественной конвекции продемонстрировано на задачах с различным значением числа Рэлея. Проведено сравнение применения различных подходов к описанию турбулентности. Анализ экспериментальных данных и результатов расчетов показывает, что качественная картина течения при значениях критерия Рэлея до 1017 моделируется верно. В случае моделирования конвективного теплообмена в кубической полости при числах Рэлея ~ 109, с экспериментальными данными совпадают, как качественные, так и количественные параметры, характеризующие процесс теплообмена. Показано, что распределение безразмерной температуры и плотности теплового потока по границам не зависит от модели турбулентности и значения критерия Рэлея, что хорошо согласуется с экспериментальными данными, однако количественные характеристики для некоторых режимов имеют существенные отличия. Для экспериментов с тепловыделяющей жидкостью (COPO и BALI) с теплоизолированной верхней границей влияние применения подходов к моделированию турбулентности не особого велико. Это может быть связано с менее интенсивным внутренним течением, поскольку теплообмен происходит только с одной границей. В рамках данных режимов результаты численного моделирования с достаточной точностью согласуются с результатами экспериментальных исследований. Для задач с двумя охлаждаемыми стенками, где присутствует интенсивная циркуляция жидкости в замкнутой области, занижается величина теплового потока на нижней границе, теплоотдача на верхней границе увеличивается. В различных источниках [3, 18, 19] также отмечается высокая погрешность по величине теплоотдачи при моделировании конвекции ТВЖ, которая возрастает при увеличении числа Рэлея. Это может быть связано с возрастанием влияния допущений, принятых при проведении расчета, поскольку набор исходных данных недостаточен для точного повторения эксперимента с помощью CFD-моделирования. К допущениям относятся отличия в задании свойств жидкости, что влияет на вычисляемые интегральные характеристики, отсутствие точных данных о характере нагрева (принимается равномерное объемное энерговыделение, хотя в эксперименте оно было неоднородным) [19], отсутствуют данные о концентрации раствора соли в воде (задаются свойства чистой воды). Отсутствие точных данных об условиях охлаждения на границах, а также данных о геометрических размерах образовавшейся ледяной корки могут вносить высокую погрешность в итоговый результат. Моделирование ледяной корки представляет серьезную проблему и не позволяет учесть шероховатость поверхности, которая может сильно влиять на коэффициент теплоотдачи. Шероховатая структура способствует дополнительной турбулезации потока и интенсификации теплоотдачи, что не учитывается при проведении численных экспериментов.
В целом, пакет программ ЛОГОС на достаточно высоком уровне воспроизводит теплообмен в ТВЖ при естественной конвекции в условиях высоких чисел Ra¿ ~ 1017. Существенную роль играет выбор подхода к моделированию турбулентности. Для режимов с охлаждаемой верхней границей, когда течение интенсивно, применение вихреразрешаю-
щей модели LES, как с качественной, так и с количественной стороны уточняет описание процессов естественной конвекции и тепломассопередачи.
Работа выполнена при финансовой поддержке гранта Президента Российской Федерации по государственной поддержке научных исследований молодых российских ученых-докторов наук МД-4874.2018.9 и ведущих научных школ Российской Федерации НШ-2685.2018.5.
Библиографический список
1. Sehgal, B.R. Core melt pressure vessel interactions during a light water reactor severe accident (MVI) / B.R. Sehgal, T.N. Dinh, R.R. Nourgaliev, V.A. Bui, J. Green, G. Kolb, A. Karbojian, S.A. Theerthan, A. Gubaidulline // FISA-97-EU research on severe accidents, EC, Luxembourg, 997. - P. 83-92.
2. Theofanous, T.G. Natural convection experiments in a hemisphere with Rayleigh numbers up to 105 / T.G. Theofanous, C. Liu // Proceedings, 1995 ANS Nat. Conf. on Heat Transfer, Portland, Oregon, 1995.
- P. 349-365.
3. Helle, M. COPO II-Lo Experiments, IVO Power Engineering LTD, YDIN-GTI-43 / M. Helle, O. Ky-malainen, E. Pessa // 4th PCRD of the European Community, 1997.
4. Kymalainen O. Heat flux distribution from a volumetrically heated pool with high Rayleigh number / O. Kymalainen, H. Tuomisto, O. Hongisto, T.G. Theofanous // Nuclear Engineering and Design. - V. 149, 1994. - P. 401-408.
5. Helle, M. Experimental data on heat flux distribution from a volumetrically heated pool with frozen boundaries / M. Helle, O. Kymalainen, H. Tuomisto // In-vessel core debris retention and coolability. Workshop proceedings, Garching, Munich, Germany, 1998. - P. 173-183.
6. Bonnet, J.M. Thermal hydraulic phenomena in corium pools: the BALI experiment / J.M. Bonnet, J.M. Seiler // Proceedings of the 7th International conference on nuclear engineering, Tokyo, Japan, 1999.
7. Bernaz, L. Thermalhydraulic phenomena in corium pools: numerical simulation with TOLBIAC and experimental validation with BALI / L. Bernaz, J.M. Bonnet, B. Spindler, C. Villermaux // In-vessel core debris retention and coolability. Workshop proceedings, Garching, Munich, Germany, 1998. - P. 185193.
8. Козелков, А.С. Исследование потенциала суперкомпьютеров для масштабируемого численного моделирования задач гидродинамики в индустриальных приложениях / А.С. Козелков, В.В. Куру-лин, С.В. Лашкин, Р.М. Шагалиев, А.В. Ялозо // Журнал вычислительной математики и математической физики. - 2016. - Т. 56. - № 8. - С. 1524-1535.
9. Kozelkov, A. Comparison of convective flux discretization schemes in detached-eddy simulation of turbulent flows on unstructured meshes / A. Kozelkov, V. Kurulin, V. Emelyanov, E. Tyatyushkina, K. Volkov // Journal of Scientific Computing. 2016. - V. 67. - P. 176-191.
10.Большухин, М.А. Об экспериментальных тестах (бенчмарках) для программных пакетов, обеспечивающих расчет теплообменников в атомной энергетике / М.А. Большухин, А.Ю. Васильев, А.В. Будников, Д.Н. Патрушев, Р.И. Романов, Д.Н. Свешников, А.Н. Сухановский, П.Г. Фрик // Вычислительная механика сплошных сред. - 2013. - Т. 5. - № 4. - С. 469-480.
11.Nourgaliev, R.R. Modeling and Analysis of Heat and Mass Transfer Processes during in-Vessel Melt Progression Stage of Light Water Reactor (LWR) Severe Accidents / R.R. Nourgaliev. - Stockholm: Royal Institute of Technology, 1998, - 458 p.
12.Fukasawa, M.Thermal-hydraulic Analysys for Inversely Statifed Molten Corium in Lower Vessel / M. Fukasawa, S. Hayakawa, M. Saito // J. Nuclear Science and Technology. - 2008. - V. 45. - № 9. -P. 873-888.
13. Abe, K. A new turbulence model for predicting fluid flow and heat transfer in separating and reattaching flows - I. Flow field calculations / K. Abe, T.Kondoh // Int. J. Heat Mass Transfer. - 1994. - V. 37. -№ 1. - P. 139-151.
14. Abe, K. A new turbulence model for predicting fluid flow and heat transfer in separating and reattaching flows - II. Thermal field calculations / K. Abe, T. Kondoh // Int. J. Heat Mass Transfer. - 1995. - V. 38.
- №. 8. - P. 1467-1481.
15.Козелков, А.С. Исследование применения RANS-моделей турбулентности для расчета неизотермических течений с низкими числами Прандтля / А.С. Козелков, А.А. Куркин, В.В. Курулин,
М.А. Легчанов, Е.С. Тятюшкина, Ю.А. Циберева // Известия РАН. Механика жидкости и газа. -2015. - № 4. - С. 44-58.
16.Menter, F.R. Ten Years of Industrial Experience with the SST Turbulence Model / F.R. Menter, M. Kuntz, R. Langtry // Turbulence, Heat and Mass Transfer 4, ed: K. Hanjalic, Y. Nagano, and M. Tummers, Begell House, Inc., 2003. - P. 625-632.
17.Menter, F.R. Zonal two-equation k-ю turbulence models for aerodynamic flows / F.R. Menter // AI-AAPaper 1993-2906.
18.Menter, F.R. Explicit algebraic Reynolds stress models for anisotropic wall-bounded flows / F.R. Menter, A.V. Garbaruk, Y. Egorov // Proc. of 3rd European Conference for Aero-Space Sciences (EUCASS). - Versailles, 2009.
19. Козелков, А.С. Моделирование турбулентных течений с использованием алгебраической модели рейнольдсовых напряжений с универсальными пристеночными функциями / А.С. Козелков, В.В. Курулин, О.Л. Пучкова, С.В. Лашкин // Вычислительная механика сплошных сред. - 2014. -Т. 7. - № 1. - С. 40-51.
20.Козелков, А.С. Зонный RANS- LES-подход на основе алгебраической модели рейнольдсовых напряжений / А.С. Козелков, О.Л. Крутякова, А.А. Куркин, В.В. Курулин, Е.С. Тятюшкина // Известия РАН. Механика жидкости и газа. - 2015. - № 5. - С. 24-33.
21.Козелков, А.С. Численная схема для моделирования турбулентных течений несжимаемой жидкости с использованием вихреразрешающих подходов / А.С. Козелков, В.В. Курулин // Журнал вычислительной математики и математической физики. - 2015. - Т. 55. - № 7. - С. 1255-1265.
22.Deryugin, Yu.N. Validation Results for the LOGOS Multifunction Software Package in Solving Problems of Aerodynamics and Gas Dynamics for the Lift-Off and Injection of Launch Vehicles / Yu.N. Deryugin, R.N. Zhuchkov, D.K. Zelenskiy, A.S. Kozelkov, A.V. Sarazov, N.F. Kudimov, Yu.M. Lipnickiy, A.V. Panasenko, A.V. Safronov // Mathematical Models and Computer Simulations. -2015. - V. 7. - №. 2. - P. 144-153.
23.Betelin, V.B. Mathematical simulation of hydrogen-oxygen combustion in rocket engines using LOGOS code / V.B. Betelin, R.M. Shagaliev, S.V. Aksenov, I.M. Belyakov, Yu.N. Deryuguin, A.S. Kozelkov, D A. Korchazhkin, V.F. Nikitin, A.V. Sarazov, D.K. Zelenskiy // Acta Astronautica. - 2014. - V. 96. -P. 53-64.
24.Kozelkov, A.S. Landslide-type tsunami modelling based on the Navier-Stokes Equations / A.S. Kozelkov, A.A. Kurkin, E.N. Pelinovsky, E.S. Tyatyushkina, V.V. Kurulin, N.V. Tarasova // Science of tsunami Hazards. - 2016. - V. 35. - №. 3. - P. 106-144.
25.Kozelkov, A.S. Numerical modeling of the 2013 meteorite entry in Lake Chebarkul, Russia / A.S. Kozelkov, A.A. Kurkin, E.N. Pelinovsky, V.V. Kurulin, E.S. Tyatyushkina // Nat. Hazards Earth Syst. Sci. - 2017. - V. 17. - P. 671-683.
26.Козелков, А.С. Численное моделирование свободного всплытия пузырька воздуха / А.С. Козелков, А.А. Куркин, В.В. Курулин, С.В. Лашкин, Н.В. Тарасова, Е.С. Тятюшкина // Известия РАН. Механика жидкости и газа. - 2016. - № 6. - С. 3-14.
27.Александров, А.А.Таблицы теплофизических свойств воды и водяного пара / А.А. Александров, Б.А. Григорьев. - М.: Издательство МЭИ, 1999, - 167 с.
Дата поступления в редакцию: 23.05.2019
S.M. Dmitriev1, O.L. Krutyakova2, A.S. Kozelkov1-3, A.A.Kurkin1, V.V. Kurulin2, D.A. Utkin2
APPLICATION OF SEMIEMPIRIC MODELS OF TURBULENCE FOR MODELING OF TURBULENT CONVECTION
Nizhny Novgorod State Technical University n.a. R.E. Alekseev1 Russian Federal Nuclear Center - All-Russian Research Institute of Experimental Physics2 Sarov Branch Institute of Physics and Technology3 of National Research Nuclear University MEPHI
Purpose: In this paper we consider the prospect of using different turbulence modeling approaches as applied to developed convection at high-order Rayleigh numbers.
Design/methodology/approach: Physical-mathematical and computational models for the simulation developed convection at height Rayleigh numbers.
Results: We demonstrate that the LES eddy-resolving model provides higher accuracy in modeling natural convection at Rayleigh numbers between 109 and 1017 than the RANS turbulence models. This advantage becomes particularly clear for the cases of Rayleigh-Bénard convection, when a large high-intensity convection zone forms at a vertical temperature difference. We demonstrate the performance of the EARSM Reynolds stress model as applied to natural convection flows in regions with dihedral corners and prevailing secondary currents. As applied to lower-intensity convec-tive flows, the differences between the involved turbulence modeling approaches are less significant. We show that the error in determining thermal hydraulic characteristics grows with increase in the Rayleigh number, and it is reasonable to use eddy-resolving turbulence modeling approaches to determine them more accurately. Area of applicability: The obtained results may be useful for the simulation of faded giant.
Key words: highly developed convective flow, turbulence modeling, high Rayleigh numbers, Rayleigh-Bénard convection, heat transfer.