CC BY
28
the impurity jet depends on the moment of momentum of the gases, the thermal power of the speed of the wind blowing the jet, and the thermophoresis conditions. The paper presents comparative graphs for calculating the level of the initial rise of the gas jet at different speeds of gas exit from the pipe and a constant wind speed.
Key words: concentration of harmful substances, emissions, gas movements, temperature, effective pipe height, working area, jet rise.
Chukarin Alexander Nikolaevich, doctor of technical sciences, professor, senior researcher, a.chukarin@ inbox.ru, Russia, Rostov-on-Don, Don State Technical University,
Manokhin Vyacheslav Yakovlevich, doctor of technical sciences, professor, ma-nohinprofa mail.ru, Russia, Voronezh, Voronezh State Technical University,
Ivanova Irina Aleksandrovna, candidate of technical sciences, docent, ivanova-eco@mail.ru, Russia, Voronezh, Voronezh State Technical University,
Golovina Elena Ivanovna, candidate of technical sciences, docent, u00111@vgasu.vrn.ru, Russia, Voronezh, Voronezh State Technical University
УДК 533.6.011.6
DOI: 10.24412/2071-6168-2022-4-263-272
ОСОБЕННОСТИ ЧИСЛЕННОГО РАСЧЕТА АЭРОДИНАМИЧЕСКОГО НАГРЕВА
Д.В. Сладков
В статье рассматривается возможность решения задачи аэродинамического нагрева осесимметричного летательного аппарата путем численного моделирования процессов газодинамического обтекания потоком воздуха и нестационарного температурного поля обтекателя. Определяется относительная точность вычислений программных комплексов Оая2, Тегт2 при решении данного рода задач.
Ключевые слова: аэродинамический нагрев, программный комплекс, сверхзвуковой поток, летательный аппарат.
В настоящее время одной из важнейших задач выступает создание отечественного программного обеспечения (ПО), позволяющего решать различные инженерные и исследовательские задачи на уровне всемирно признанных вычислительных комплексов. Особое значение при этом имеют CAE-системы, дающие возможность методами численного моделирования осуществлять симуляцию течения физических процессов и таким образом заменять часть затратных натурных экспериментов, необходимых при разработке новых изделий.
При разработке такого ПО крайне важно определение точности вычислений путем сравнения получаемых результатов с решениями, полученными при использовании хорошо зарекомендовавших себя аналогах и экспериментальными данными. Важной особенностью выступает проведение таких исследований для каждого рода решаемых задач, поскольку математическое описание моделируемых процессов может значительно отличаться в зависимости от рассматриваемых условий, геометрии расчетной области и объекта исследований, применяемых веществ и материалов.
Одним из актуальных направлений является решение задач аэродинамического нагрева летательных аппаратов (ЛА) при полете с высокими сверхзвуковыми скоростями. Необходимо отметить, что температура ЛА экспериментально определяется крайне редко. Это связано с высокой стоимостью и сложностью проведения такого рода исследований, поэтому гораздо чаще используются аналитические вычисления или численное моделирование.
В настоящей работе применяется следующий подход для решения задачи аэродинамического нагрева [1]:
1. Проводится расчет траектории полета ЛА и определение профиля скорости. На основе полученных данных выбираются ключевые точки, для каждой из которых определяются значения скорости полета, высоты, параметры атмосферы.
2. Решается газодинамическая задача обтекания ЛА в ключевых точках траектории, по результатам которой формируется набор данных, составляющих граничные условия тепловой задачи: функции локального коэффициента теплопередачи и полной температуры газа по длине стенки ЛА.
3. Осуществляется расчет нестационарного температурного поля элементов планера ЛА на основе полученных параметров газа вблизи стенки корпуса для ключевых точек траектории. Промежуточные значения тепловых характеристик газа определяются методом линейной интерполяции.
Исследование проводится на базе программных комплексов (ПК) Gas2 (газодинамика) и Term2 (тепловая задача), разработанных в ТулГУ, для осесимметричного ЛА, а точность вычислений определяется относительно ПК Ansys, в котором для решения газодинамической и тепловой задач применяются модули Fluent и Thermal соответственно.
Первоначально необходимо верифицировать модель обтекания ЛА потоком воздуха. Для этого зачастую применяется сравнение рассчитанных в различных комплексах аэродинамических характеристик (АДХ) ЛА, одной из которых является безразмерный коэффициент лобового сопротивления Cx, зависящий в основном от скорости полета и определяемый по зависимости [2]:
Cx Q ■S,
где F- сила, воздействующая на ЛА, Q — р-v2 /2 - скоростной напор, S- площадь миделя.
В этой связи требуется провести серию расчетов в рассматриваемых ПК по определению зависимости Cx от скорости полета. Важно отметить, что такие задачи зачастую решаются в обращенной постановке - расчетную схему строят таким образом, чтобы воспроизвести процесс продувки в аэродинамической трубе, когда объект исследования остается неподвижным, а поток воздуха обтекает модель со скоростью, равной скорости полета.
Вычисления в ПК Ansys Fluent и Gas2 осуществляются с применением к—£ модели турбулентности, включающей использование стандартных пристеночных функций [4]. При этом важной особенностью Gas2 является определение силы, воздействующей на тело, лишь на основе давления - без учета силы вязкого трения, что необходимо учитывают при расчетах.
В качестве объекта исследования для решения тестовых задач принята схема осесимметричного ЛА с диаметром фюзеляжа 0,29 м и длиной 3,7 м. Расчеты проведены при параметрах околоземного атмосферного давления, равного 101325 Па, и температуре окружающей среды 293 К. Набор скоростей набегающего потока воздуха определен исходя из значения местной скорости звука и соответствует следующим числам Маха: 1; 1,5; 2; 2,5; 3; 3,5; 4 М
Сетка дискретизации расчетной области с приближением в область носовой части ЛА представлена на рис. 1.
[ШШЩШШПШШШ
У Г 1.4 Г gj Г Г Г-1 -IГ-1 jrrl.rrl.j-J Г"! jrrljrijrrrl.rrl.-ir Г-1Г Г Г"! -Д Г-1 jrrl.-irl-l'rr г4 Г-1 jri.-irri.rrl.-I-i Г"! jrrijfijsrrl.rfl.-is fj Г Г Г-1 -IГ-1 jrrl.-irl.j4r г-д rrljrijrrrirrl.-irrl.rr Г"! -Д Г"! jrrijrl-lzrrl.rfl.-ir Г-Д Г Г Г-1 j r-1 "j г г
ЩщШ|:11: i i ri j: шг liimmhrnfптШНпгптШ^ j | И И ff
llllllllllllllllllll JIM
ШгНгшЗпШШ
а
б
Рис. 1. Сетка дискретизации расчетной области: а - созданная в Vis2, б — созданная в Ansys Mesh
На рис. 2 приведены зависимости Cx от числа Маха, рассчитанные без учета силы вязкого трения, воздействующей на ЛА, полученные в Gas2 и Ansys Fluent.
HP4
Cx 0,45 0.4 0,35 0,3 0,25 0.2 0,15 0,1
1 1,5 2 2,5 3 3,5 M
Рис. 2. Зависимость коэффициента лобового сопротивления Cx ЛА от числа Маха M: 1 — результаты расчета в ANSYS Fluent, 2 - результаты расчета в Gas2
Отклонения результатов, полученных в программном комплексе Gas2, от величин, полученных в комплексе Ansys Fluent, представлены в табл. 1.
Относительные ошибки расчетных (ANSYS Fluent и Gas2) Cx
Таблица 1
M 1 1,5 2 2,5 3 3,5 4
А, % -8,69565 0 -2,85714 -12,9032 -14,8148 -12 -13,0435
Из рис. 2 и табл. 1 видно, что ПК Gas2 дает заниженные значения Сх, при этом максимальные отклонения наблюдаются при М > 2,5 и достигают 14,8%, что неизбежно будет сказываться на точности моделирования нагрева корпуса ЛА.
Важно отметить, что при решении такого рода задач значительное влияние на точность результатов оказывает степень и вид дискретизации расчетной области в зоне пограничного слоя ЛА. Поскольку в программе Vis2 используется метод дискретизации на ячейки одинакового размера во всей области, точно разрешить пограничный слой не представляется возможным. Это оказывает существенное влияние на определение параметров газа в указанной области, в особенности при высоких сверхзвуковых скоростях и в случае моделирования сложных процессов - переход из дозвуковой в сверхзвуковую скорость.
На рис. 3, а представлена траектория полета ЛА, рассчитанная для решения тестовых задач, при горизонтальном угле старта 75°, а на рис. 3, б показан закон изменения скорости полета во времени - профиль скорости, полученные в разработанной и зарегистрированной программе [6]. Параметры атмосферы для различных высот принимались в соответствии с ГОСТ 4404-81.
Н Ье^тркнрмпмнКи-КТржчдоТчм».!)^**" -ПК
ОДкПЫ НКЧ±Т СШ ГС№Ш
о е н к х | я Ч & I # 9 и з
^оннив Грвфцхн результате« 30 Грвфики результат» 20 Текст рвзугынтов Настрой*«
ито» ¡ПК 4«1»ЙИ зга*'« г*а*'0» 15*еИИ •СЪечО* ££5М>0* Хм
тм ЛЧ УрО Дмм ийпе V Ю " х Ю * а 0 Х- 700000 « «н « 0 -«у-с 100000
и .ригп^нн ГНЛНЧ РС ¡4 ■ СЛТ^фт^рЛТв^имОИ&.ТмрЛл - О
ФнЯлы Р*счст Еид ГТоиощ»
□ вн ► | х яяя^ез х| ? из
Донмие Греф»« результатов 30 Графики р09ульт«тф& 20 Текст роэультвточ Ностройад
Рис. 3. Типовой закон изменения скорости полета во времени: а — траектория полета ЛА; б — профиль скорости
266
Угол пуска был выбран для проведения серии расчетов при различных начальных и граничных условиях, соответствующих как плотным (вблизи поверхности земли), так и разреженным слоям (стратосфера) атмосферы. Это способствует более полной отработке связки программных комплексов Gas2-Term2.
Расчет газодинамического обтекания проводится для следующего набора точек траектории:
- высота полета 0 м, температура воздуха 293 К, время t=0 с;
- промежуточные значения между стартовой и максимальной скоростями, высота полета 1000, 1350 и 2600 м, температура воздуха 293, 279 и 271 К, время 2, 3,7, 5 с;
- высота полета 3650 м, температура воздуха 264 К, время 5,8 с;
- высота полета 12800 м, температура воздуха 216 К, время 13,2 с;
- высота полета 63860 м, температура воздуха 237 К, время 111,8 с;
- высота полета 15000 м, температура воздуха 216 К, время 215 с;
- высота полета 0 м, температура воздуха 293 К, время 232,9 с.
На рис. 4 приведены примеры картин распределения температур в расчетной области для максимальной скорости полета ЛА - 1540 м/с. При этом для результатов, полученных в ANSYS Fluent диапазон палитры приведен к значениям, полученным в Gas2, что необходимо для более наглядного представления результатов.
Рис. 4. Поля температуры, полученные в Gas2 (верхний рисунок) и в Ansys Fluent и соответствующие точке траектории в момент времени 5,8 с
Из приведенных рисунков видно, что полученные в ПК распределения температур отличаются незначительно, отрывные течения со скачками уплотнения близки по характеру и углу отрыва.
Основные различия наблюдаются в зоне пограничного слоя, что обусловлено несколькими факторами:
1. Препроцессор ПК ANSYS при построении сетки позволяет проводить локальное измельчение, за счет чего происходит разрешение пограничного слоя.
267
2. В качестве граничного условия на стенке корпуса ЛА в Gas2 задается условие полного проскальзывания частиц, в ANSYS Fluent задается условие торможения, в связи с чем температура пограничного слоя во втором случае равна температуре торможения потока.
Моделирование газодинамического процесса обтекания корпуса позволяет определить параметры теплообмена на его поверхности, выступающие граничными условиями тепловой задачи при рассчете температурного поля обечайки. Последовательная серия расчетов на множестве точек траектории полета позволяет сформировать набор граничных условий, импортируемых в ПК Term2 для определения нестационарного температурного поля корпуса ЛА.
В качестве материала обечайки принята сталь 30 без теплозащитного покрытия с толщиной стенки 2,2 мм.
На рис. 5 показан контур исследуемой на тепловое состояние обечайки, являющейся расчетной областью, и границы теплообмена для расчета аэродинамического нагрева (красным цветом показана граница, нагреваемая обтекающим потоком, зеленым цветом - граница контакта с воздухом внутренней полости).
D а р a IbQLQtt ЪВ1Н9В|!№ ЧЧЧ Щ f
1 1 Mel @ ^ llfИЛА Ф Q Ш
□ н v 1 ш ш l В » и ш Ш ч. & ..,,, ?
Дискретизация данной области на конечные элементы, а также точки для построения графиков изменения температуры во времени, приведены на рис. 6. Необходимо отметить, что приведенные контур расчетной области, границы теплообмена, а также расположение точек снятия температуры для ПК Ansys аналогичны.
Результаты расчетов нестационарного температурного поля корпуса в выбранных точках приведены на рис. 7 - 8.
Характер кривых изменения температуры во времени в значительной степени определяется не только законом изменения скорости, но и изменением высоты полета и плотности воздуха на траектории. Наиболее нагретой является носовая часть обтекателя, что вызвано сильным торможением потока в данной зоне. Температура на носке достигает 900 К и более. Максимальные температуры возникают на 8-10 секундах полета, когда ЛА движется в плотных слоях атмосферы со скоростью близкой к максимальной.
Также на рис. 8 видно, что характер кривых, полученных в Term2, не имеет значительных отличий от результатов из Ansys Thermal, что подтверждает правильность не только математической модели, заложенной в программный комплекс, но и распределения тепловых параметров на границе.
На основании полученных кривых изменения температуры был проведен сравнительный анализ значений температур в каждой точке для ключевых моментов времени. Полученные значения относительной погрешности вычислений по программе Term2 относительно значений, полученных с помощью комплекса Ansys Thermal сведены в табл. 2.
В ш|~в ■ « ft, % Ч J f
iffitifeiiiii ^вд-^тГЩ^ as idtaai' sis Ф oi gj
Параметры: n4u=63551 nu =57557
IdU =139
IdE =136 Тот
Рис. 6. Расчетная сетка, созданная в Vis2 (Верхнее изображение) и в Ansys Mesh, расположение точек корпуса, для которых выполняется построение графиков изменения во времени температуры стенки
б
Рис. 7. Поле температуры стенки, полученное для моментов времени t=8,3 c и t=232,9 c, где а — результаты вычислений в Term2, б — в Ansys Thermal
269
Рис. 8. Графики изменения температуры стенки, полученные в Term2 (верхний рисунок) и в Ansys Thermal
Таблица 2
Относительные ошибки расчетных температур по результатам вычислений
Относительные ошибки, %
№ точки 5,8 8,3 232,86
1 -10,9521 -7,23581 -9,33582
2 -12,0397 -12,0552 -8,69223
3 -11,2264 -8,3007 -1,30756
4 -11,1976 -10,3089 -5,50189
5 -10,4503 -11,1773 -6,19797
6 -5,32731 -4,81492 -3,83021
7 -1,90015 -4,00221 -3,52772
8 -10,4791 -10,8815 -7,23336
Из приведенных данных видно, что показатель относительных ошибок находится в диапазоне от 1,31 до 10,95 %, причем для большинства точек корпуса наибольшие отклонения возникают в момент времени, соответствующий максимальной скорости полета ЛА. Результаты расчета температуры в Term2 имеют заниженное значение, что связано с невозможностью точно разрешить пограничный слой при равномерной дискретизации расчетной области. Это может приводить к критичным ошибкам при проектировании, поскольку закладываемый запас прочности и устойчивости к тепловым нагрузкам у изделия может оказаться недостаточным.
Следует отметить, что аналогичная ситуация наблюдалась при верификации газодинамической модели обтекания (рис. 2, табл. 1), когда относительные отклонения были наибольшими при высоких сверхзвуковых скоростях. Полученная разница в результатах расчетов температуры находится в тех же пределах, что ошибки, полученные при расчете АДХ изделия, что свидетельствует о правильности математической модели нестационарного нагрева, а наиболее вероятной причиной возникновения ошибок в тепловом расчете являются неточности в определении аэродинамических параметров при решении газодинамической задачи.
При этом в большинстве случаев относительные погрешности в вычислениях не превышают 10%, что является приемлемым результатом, в особенности - при решении инженерных задач, когда помимо точности решения большое значение имеет время, затрачиваемое на вычисления.
В дальнейшем необходима разработка новых методов дискретизации расчетной области, которые должны включать возможность локального измельчения сетки, поскольку это даст возможность разрешать пограничный слой ЛА, а также учесть при расчетах влияние силы вязкого трения. Как следствие, это позволит исключить необходимость в разработке методик, основанных на применении критериальных зависимостей и перейти к сопряжению газодинамических и тепловых задач посредством граничных условий четвертого рода.
Список литературы
1. Авдуевский В.С., Данилов Ю.И., Кошкин В.К. и др. Основы теплопередачи в авиационной и ракетно-космической технике. М.: Машиностроение, 1975. 624 с.
2. Дмитриевский А.А., Лысенко Л.Н. Внешняя баллистика: учебник. М.: Машиностроение, 2005. 608 с.
3. Лыков А.В., Алексашенко А.А., Алексашенко В.А. Сопряженные задачи конвективного теплообмена. Минск: Изд. Белорус. ун-та, 1971. 346 с.
4. Новожилов В.В. Теория плоского турбулентного пограничного слоя несжимаемой жидкости. Л.: Изд. Судостроение, 1977. 165 с.
5. Патанкар С. Численные методы решения задач теплообмена и динамики жидкости. М.: Энергоиздат, 1984. 152 с.
6. Свидетельство 2021617512. Программа расчёта начальных участков движения РС оснащённых стартовым двигателем: программа для ЭВМ / В.А. Дунаев, Д.Б. Скорлупкин, В.Ю. Сладков, Д.В. Сладков (RU); правообладатель Акционерное общество «Научно-производственное объединение «СПЛАВ» имени А. Н. Ганичева» (АО «НПО «СПЛАВ» им. А. Н. Ганичева»). № 2021614909; заявл. 06.04.21; опубл. 17.05.2021, Бюл. № 5, 18,5 Мб.
Сладков Дмитрий Валерьевич, магистрант, sladckov. d@yandex. ru, Россия, Тула, Тульский Государственный Университет,
Научный руководитель - Дунаев Валерий Александрович, д-р техн. наук,
профессор, dwa222@mail. ru, Россия, Тула, Тульский государственный университет
FEATURES OF NUMERICAL CALCULATION OF AERODYNAMIC HEATING
D.V. Sladkov
The article considers the possibility of solving the problem of aerodynamic heating of an axisymmetric aircraft by numerical simulation of the processes of gas-dynamic flow around the airflow and the non-stationary temperature field of the fairing. The relative accuracy of calculations of the software systems Gas2, Term2 is determined when solving this kind ofproblems.
Key words: aerodynamic heating, software package, supersonic flow, aircraft.
Sladkov Dmitri Valeryevich, student, sladckov.d@yandex.ru, Russia, Tula, Tula State University,
Scientific adviser - Dunaev Valery Alexandrovich, doctor of technical sciences, professor, dwa222@mail.ru, Russia, Tula, Tula State University
УДК 004.94:502.1
DOI: 10.24412/2071-6168-2022-4-272-275
МЕТОД ИССЛЕДОВАНИЯ РАБОТЫ ГРАДИРНИ КОМПЬЮТЕРНЫМ МОДЕЛИРОВАНИЕМ И АНАЛИЗОМ ДАННЫХ
А.М. Бочарова
В статье затрагивается тема о том, как влияет величина элемента при расчете в программных комплексах, основанных на методе конечных элементов на некоторые параметры, в числе которых качество расчета и затраченное на моделирование время. Дается сравнение компьютерного моделирования работы градирни при различных условиях. Рассматривается несколько вариантов разбиения сетки конечных элементов, а именно в зависимости от размера одного элемента. Особое внимание уделено оценке влияния количества и размеров отдельных элементов на различные параметры процесса, такие как скорость выдуваемого из градирни воздуха, его температура, смешение воздуха из градирни с окружающей средой, и нескольких других параметров. Целью статьи является сравнительный анализ параметров нескольких вариантов компьютерных моделирований, которые отличаются между собой исключительно сеткой конечных элементов. Приводятся конкретные значения оцениваемых величин в зависимости от размера одного элемента. Автор приходит к выводу, о том, какой вариант наиболее предпочтительный и оптимальный с точки зрения погрешности расчета и его скорости.
Ключевые слова: анализ, информация, данные, программное обеспечение, программа, компьютерное моделирование, экология.
Разнообразное программное обеспечение закрепилось в жизни современной науки, она помогает во многих видах научной деятельности: автоматизация расчетов по известным формулам, оцифровка карт, подготовка рукописей, коммуникация с коллегами, проверка орфографии и заимствований, патентный поиск, исследование литературы, а также большое значение имеют программные комплексы для решения прикладных физических, химических и иных задач [1-8].
272
