Физические особенности обтекания датчиков обледенения и решение практических задач, связанных с полетом гражданского самолета в условиях обледенения тема диссертации и автореферата по ВАК РФ 05.07.01, кандидат наук Нгуен Нгок Шанг

• Нгуен Нгок Шанг
• кандидат наук
• 2022

Оглавление диссертации кандидат наук Нгуен Нгок Шанг

Введение

Глава 1. Обзор литературы

1.1 Определение аэродинамических параметров обтекания ЛА

1.1.1 LEWICE и ONERA

1.1.2 FENSAP-ICE

1.2 Моделирование движения переохлажденных капель

1.2.1 LEWICE и ONERA

1.2.2 FENSAP-ICE

1.3 Моделирование образования льда

1.3.1 Модель Мессингера

1.3.2 Мелководная модель SWIM

1.4 Итог обзора литературы

Глава 2. Математическая постановка и численный метод

2.1 Метод конечного объёма

2.1.1 Шаблон численной схемы

2.1.2 Градиенты в центрах ячеек

2.1.3 Ограничители градиентов в центрах ячеек

2.1.4 Градиенты на гранях ячеек

2.2 Моделирование течения воздуха

2.2.1 Эффект шероховатости

2.2.2 Модель Спаларта-Аллмараса с учетом шероховатости

2.2.3 Применение закона стенки

2.3 Модели движения переохлажденных капель

2.3.1 Системы уравнений

2.3.2 Метод расчета потоков через грань ячейки

2.3.3 Граничное условие на поверхности тела

2.3.4 Распределение капель в облаке

2.3.5 Моделирование движения ледяных кристаллов

2.4 Модель термодинамического процесса нарастания льда

2.4.1 Численный метод решения системы SWIM

2.4.2 Нарастание льда

Глава 3. Верификация разработанной программы

3.1 Солвер водности

3.1.1 Одномерное движение капель

3.1.2 Двумерное движение капель

3.1.3 Двумерное движение ледяных кристаллов

3.2 Модификация солвера воздуха

3.2.1 Эксперименты Хосни

3.2.2 Эксперименты Хеалзера

3.3 Солвер термодинамического процесса нарастания льда

3.3.1 Обледенение профиля NACA0012 (Glaze Ice)

3.3.2 Обледенение профиля NACA0012 (Rime Ice)

Глава 4. Физические особенности обтекания датчика обледенения

4.1 Принцип работы вибрационного сигнализатора обледенения

4.2 Нарастание льда на датчике обледенения в неравномерном поле водности

4.3 Выбор места расположения датчиков обледенения

4.4 Методология определения исходных данных для стендовых испытаний датчиков системы воздушных сигналов в условиях обледенения

Заключение

Список литературы

Введение диссертации (часть автореферата) на тему «Физические особенности обтекания датчиков обледенения и решение практических задач, связанных с полетом гражданского самолета в условиях обледенения»

Введение

Обледенение - это процесс образования льда на поверхностях различных предметов при низких температурах. Это явление сильно и, зачастую, негативно влияет на условия жизни и деятельности человека, особенно в тех местах, где зима суровая. При температурах ниже нуля вода и водяной пар замерзают и превращаются в лёд. Вследствие мороза, возможны отказы в работе агрегатов, неустойчивость или неисправность машин, что может приводить к авариям.

С первых же полётов летательных аппаратов (ЛА), люди заметили влияние обледенения на них. В полете, когда ЛА летает в облаках, содержащих переохлажденные (температура ниже нуля) водяные капли, может происходить обледенение обтекаемой поверхности после столкновения капель с поверхностью. Вследствие этого происходит обледенение деталей самолета: крыльев, передней части фюзеляжа, воздухозаборников, воздушных винтов, элементов оперения. Ледяной слой изменяет форму обтекаемых поверхностей, портит их аэродинамическую форму и, как следствие, уменьшает подъемную силу, увеличивает сопротивление, ухудшает управляемость и устойчивость. Кроме этого, обледенение может вызвать неточность или даже отказ в работе датчиков, приборов, агрегатов, и нарушить режим работы двигателя.

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

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

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

Таким образом, исследование обледенения является важной и актуальной задачей. Для изучения этого процесса, используют три метода: теоретический, экспериментальный и компьютерное моделирование. Если два первых метода применяются достаточно давно, то третий метод появился только с начала восьмидесятых годов. Эти методы не существуют независимо, они дополняют друг от друга и помогают друг другу совершенствоваться.

Обледенение представляет собой очень сложное физическое явление. Оно состоит из ряда разных, но зависимых друг от друга процессов, таких как обтекание летательного аппарата воздухом, движение водяных капель,

взаимодействие капель с поверхностью тела, теплопередача, образование льда и т.д. Поэтому невозможно создать универсальную теорию для всех случаев обледенения. Используя теоретический метод, мы можем решать ограниченное число простых задач обледенения. У экспериментального метода возможности также невелики, поскольку нужно создать реальные условия обледенения и точную модель объекта, либо нужно использовать оригинальный объект. Это очень сложно и стоит дорого! Поэтому метод компьютерного моделирования благодаря своим уникальным возможностям стал лидером в изучении обледенения. Это прежде всего гибкость в изменении внешних условий обледенения и формы объекта. С развитием компьютерной науки, компьютерное моделирование становится важной и необходимой частью в любом научном исследовании, в том числе исследовании обледенения.

В настоящее времени существует много коммерческих пакетов прикладных программ вычислительной гидродинамики, таких как FlowVision [1], ANSYS [2], NUMECA [3], COMSOL Multiphysics [4], Autodesk CFD [5],... и программ с открытым исходным кодом, таких как OpenFoam [6], SU2 [7]. Несмотря на это, только некоторые программы поддерживают моделирование обледенения, такие как FlowVision [8], ANSYS FENSAP-ICE [9]. В этом же ряду находятся исследовательские коды LEWICE [10], ONERA3D [11], NSCODE-ICE [12], EWT-ЦАГИ [13].

Целью данной работы являются:

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

Для достижения поставленной цели решены следующие задачи:

1. Проведен аналитический обзор источников литературы и выбран метод расчета процесса нарастания льда на поверхностях самолета с учетом шероховатости поверхности и мелкодисперсного состояния водной фракции воздушного потока.

2. Разработана эффективная методика, позволяющая надежно решать систему уравнений водности, описывающую в трехмерном приближении движение переохлажденных капель водной фракции потока.

3. Реализован метод SWIM нарастания льда на поверхности датчика обледенения и других поверхностях летательного аппарата.

4. Выполнена серия верификационных и валидационных расчетов.

5. На основе реализованной автором численной методики выполнено исследование физических особенностей обтекания датчиков обледенения.

6. Разработана методология определения исходных данных для стендовых испытаний датчиков системы воздушных сигналов в условиях обледенения на основе численного моделирования.

7. Разработана упрощенная методология предварительного выбора места расположения датчиков обледенения на основе численного моделирования физических особенностей нарастания льда в неравномерном поле водной фракции потока.

8. Решена практическая задача, связанная с полетом гражданского самолета в условиях обледенения: осуществлен выбор места расположения датчиков обледенения на поверхности самолета.

Научая новизна:

1. Впервые расчетным образом исследованы физические особенности обтекания датчика обледенения в условиях компоновки с элементами планера и сформулированы принципы выбора места расположения датчика

обледенения на поверхности самолета. В частности, показано, что скорость нарастания льда в неравномерном поле водной среды с точностью не хуже 3% совпадает со скоростью нарастания льда в «эквивалентном» равномерном поле водной среды.

2. Применение модифицированной процедуры «распада произвольного разрыва» в совокупности с использованием интегрирования методом «локального шага» по времени позволяет не менее, чем на 10% ускорить расчет полей водности по сравнению с существующими методами.

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

Методология и методы исследования. В диссертации применяются численные методы для моделирования сжимаемых течений, движения переохлажденных капель и термодинамического процесса нарастания льда. Для реализации программы используются языки программирования С++ и Python.

На защиту выносится:

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

2. Методология определения исходных данных для стендовых испытаний датчиков системы воздушных сигналов в условиях обледенения на основе численного моделирования.

3. Упрощенная методология предварительного выбора места расположения датчиков обледенения на основе численного моделирования физических особенностей нарастания льда в неравномерном поле водной фракции потока и результаты исследований физических особенностей обтекания подобных датчиков.

4. Примеры решения практических задач, связанных с полетом гражданского самолета в условиях обледенения.

Достоверность полученных результатов подтверждается сравнением с аналитическими решениями, с экспериментальными данными и с результатами других авторов.

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

Апробация работы. Основные результаты работы докладывались автором на II Международной конференции «Математическое моделирование», Москва, МАИ, 21-22 июля 2021г; 28-й Всероссийской конференции с международным участием «Высокоэнергетические процессы в механике сплошной среды», Новосибирск, ИТПМ СО РАН, 20-24 сентября 2021г; 64-й Всероссийской научной конференции МФТИ, Жуковский, МФТИ, 03 декабря 2021г; «Отечественные CFD коды — 2021», Москва, ИПМ РАН, 18-19 декабря 2021.

Публикации. Основные результаты по теме диссертации изложены в 4 научных публикациях, которые опубликованы в научном издании, включенном в списке RSCI. Список публикаций:

1. Нгуен Нгок Шанг, Численный метод для модели Эйлера движения переохлажденных капель в условиях обледенения, Труды МФТИ Том 13, № 3 (51) 2021, стр. 133-143

2. Нгуен Нгок Шанг, Граничное условие «закон стенки» с учетом шероховатости для модели турбулентности Спаларта-Аллмараса, Труды МФТИ Том 13, № 4 (52) 2021, стр. 103-113

3. Нгуен Нгок Шанг, Выбор места расположения датчика обледенения и датчика полного давления в условиях обледенения, Труды МФТИ Том 14, № 1 (53) 2022, стр. 3-14

4. Нгуен Нгок Шанг, Методология определения исходных данных для стендовых испытаний датчиков системы воздушных сигналов в условиях обледенения, Труды МФТИ Том 14, № 1 (53) 2022, стр. 15-26

Структура и объём диссертации. Данная работа состоит из введения, 4 глав, заключения и списка литературы, содержащего 63 наименований. Объём текста составляет 126 страниц, включая 84 рисунков и 14 таблиц.

Глава 1. Обзор литературы

Компьютерные программы моделирования обледенения разрабатываются во многих странах, таких как США, Франция, Канада, Англия, Италия, Россия... Первые программы, например LEWICE (США) [10, 14], О№КЛ (Франция) [11, 15], начинали разрабатывать в начале 1980 годов. В данный момент, одной из самых успешных программ является FENSAP-ICE (США-Канада) [9, 16], которая была разработана в конце 80х годов. В России в 2020г. была реализована новая версия пакета FlowVision [8], в котором появилась возможность моделирования процесса обледенения поверхностей - как на сухом, так и влажном режимах. Поэтому, главными источниками литературы являются документации этих программ и статьи их авторов.

Эта литература показала, что для моделирования обледенения нужно решать следующие задачи:

1. Моделирование течения воздуха,

2. Моделирование движения переохлажденных капель,

3. Моделирование образования льда,

4. По мере нарастания льда форма летательного аппарата меняется, приходится перестраивать расчетную сетку и циклически повторять решение первых трех задач.

Рисунок 1.1 - Связь между модулями

Связь между решениями этих задач показана на рисунке 1.1. В данной работе для построения расчетной сетки используются программы с открытым исходным кодом Gmsh [17] и Salome [18]. Поэтому обзор литературы проведен только для трех первых задач.

1.1 Определение аэродинамических параметров обтекания ЛА

Первым шагом при моделировании процесса обледенения является определение аэродинамических характеристик потока воздуха. Полученные результаты будут использоваться как входные параметры для моделирования потока переохлажденных капель и термодинамического процесса нарастания льда. Для расчета поля водности необходимы скорость и плотность воздуха; для термодинамического расчета - температура, напряжение трения на стенке и коэффициент теплообмена на поверхности.

1.1.1 LEWICE и ONERA

Поле воздуха вычисляется разными методами в разных программах. Программа LEWICE [10, 19] использует панельный метод для вычисления потенциальных течений. В ONERA [11, 15] кроме этого подхода, применяется метод конечного объёма для решения уравнений Эйлера. Эти программы моделируют невязкое течение, полученные результаты используются для моделирования поля водности. Для термодинамического моделирования обледенения нужны параметры вязкого течения, поэтому в этих программах решается уравнение пограничного слоя. Надо обратить внимание на то, что когда лёд появляется, шероховатость поверхности увеличивает коэффициент теплообмена и напряжение трения на стенке. Первой трудностью при изучении шероховатости является многообразие форм неровностей. Для решения этой проблемы была придумана некая идеальная модельная шероховатость: «эквивалентная песочная модельная шероховатость» с одинаковой высотой бугорков ks.

Rek = Refc (1)

В LEWICE требуется определять момент переход от ламинарного к турбулентному режиму течения в пограничном слое. Критерием перехода является шероховатое число Рейнольдса:

r i ~ 600 ламинарный v ' к {> 600 турбулентный где uk - скорость на высоте элемента неровности ks, которая вычисляется по

формуле:

ик 3 4 IS2 due ks

— = 2ksd — 2ksd + ksd + —— -j^-ksd(l — ksd) , ksd =

где s - расстояние до точки торможения, ue - скорость на внешней границе пограничного слоя. Толщину пограничного слоя можно определить с помощью формулы:

S « 8.5вь

где ламинарная толщина потери импульса вычисляется по формуле Туэйтса:

„ 0.45va fs _ в22 =-— I и5е ds.

Ug JQ

Тогда ламинарный коэффициент теплоотдачи можно определить, используя формулу, разработанную Смитом и Сполдингом:

1-1/2

Л

hc = 0.296-= Vv

щ288 I u188ds

'о

(2)

[ иУ ¿0

где X - теплопроводность воздуха.

Для турбулентной области коэффициент теплопередачи вычисляется по формуле:

кс = Б1риеср, (3)

где - число Стантона:

сг/2

St =

Рт

Prt+ Stk

Здесь:

Prt - турбулентное число Прандтля: Prt = 0.9,

iu к \

Stk - шероховатое число Стантона: Stk = 1.16 (

ит - динамическая скорость: ит = ие ¡у.

0.2

Cf - коэффициент трения: — =

0.41

1п(864.0^)+2.568

Турбулентная толщина потери импульса:

[0.0156 Г5

dt(s) =

и

4.11

f и

stran

ul86 ds

0.8

+ 9l(strau).

В ОМЕЯЛ уравнения пограничного слоя моделируются с использованием метода конечного объёма. Переход в турбулентный режим определяется так же как в LEWICE. Для модели шероховатости, в отличии от LEWICE, используется формула Ван Дриста. Возрастание коэффициента теплопередачи обеспечивается увеличением турбулентного напряжения трения гладкой поверхности:

xt = pF2l^

22 /ди\ /dw\

vdyj +(ду/

где длина смешения:

( к У\

1т = 8к^апк(— §)'Кг = 0.085, к = 0.41.

Турбулентный фактор демпфирования возле гладкой стенки вычисляется по формуле:

4mjrp\ ехР I ~26к/7~)'

Эффект шероховатости моделируется модификацией F в модели Ван Дриста:

(f60у\

Тогда коэффициент тепловой передачи вычисляется по формуле:

/ дТ \

К = Л ( —- J (7wall - Tree У1'

У=0

где TWaii - температура поверхности, Tree - равновесная температура.

2

1.1.2 FENSAP-ICE

FENSAP-ICE [9, 16] решает уравнения RANS методом конечного элемента. Чтобы учесть эффект шероховатости используется модификация модели турбулентности Спаларта-Аллмараса [20]. Расчеты проводятся с граничным условием постоянной температуры на стенке. По результатам расчетов определяются поток тепла на стенке и коэффициент теплопередачи:

, Qwall

hc =--——. (5)

* wall * сю

Таким образом, в отличии от LEWICE и ONERA, FENSAP-ICE использует только один солвер, один метод чтобы получить все необходимые параметры течения воздуха. Такой метод более универсален, его можно использовать для сложных геометрий и разных режимов течения.

1.2 Моделирование движения переохлажденных капель

При движении переохлажденных капель, может происходить много физических процессов, например, столкновение, разделение, соединение капель. Учет всех этих процессов приводит к существенному усложнению моделирования движения капель. Поэтому обычно делаются следующие фундаментальные допущения:

• отсутствует влияние капель на поле воздуха,

• отсутствуют столкновения, разделение или соединение капель,

• переохлажденные капли имеют форму твердой сферы.

Важными параметрами являются:

• й - диаметра капли

• ЫУО - средний диаметр капель,

• ЬЖС - водность облаков.

Размер капли будет оказывать влияние на силу инерции и силу сопротивления. Маленькие капли до попадания на поверхность ЛА будут, в

основном, двигаться по линиям тока воздуха, в то время как большие капли до столкновения с поверхностью ЛА будут, в основном, двигаться по инерции. В облаке содержатся капли разных размеров, поэтому важным фактором является средний размер капель и их распределение по размеру.

Для моделирования движения переохлажденных капель, применяются разные методы. Однако, в результате должны быть получены скорость капель до столкновения (Vd), массовой поток попадающих капель (mimp) и коэффициент захвата (ß).

1.2.1 LEWICE и ONERA

В обеих программах используется подход Лагранжа к описанию движения капель. Уравнения движения капель имеют следующий вид [10, 15]:

dXr!

—~ = Vd

dVd ^ (6)

<та^[Г = 0-5CdpaAd\\Vr\\Vr + mdg,

4 "2

где xd - положение капли, Vd - скорость капли, md = -pwnrd - масса капли, Cd -

коэффициент сопротивления, pa - плотность воздуха, Ad = nrd - характерная площадь капли, Vr = Va(xd)-Vd - скорость капли относительно скорости воздуха, д - ускорение свободного падения. В ONERA [15] коэффициент сопротивления определяется по формуле Шиллера и Наумана:

24

CD=—(l + 0.l5ReT7),

Rew

где

В LEWICE [21]:

2rd\\Vr\\ Re d = d" r" < 800. d v

24 6

CD =-+ 0.4 +-г^,если Cd < 100,

Red 1 + VRed

24 6 CD =-+ 0.3 +-г^.если Cd > 100.

Red l + VRed

Для решения лагранжевой системы уравнений движения капель используются разные вычислительные методы. В ОМЕЯЛ траектория капель вычисляется по следующей схеме:

хп+1 =хп + Уа(хп)Л1- + [(уй)п - уа(хп)}

1 — ехр

- Л

п

п

+

—ле

М — [1 — ехрх ^

а

п

49'

(7)

(Уа)П = Уа(хП) + ехр

— Л

[(Уа)П — Уа(хП)] +

1 — х р

— Л

тпа9,

24 2рыг2

где та =---. Это схема первого порядка без условия сходимости. Можно

использовать схемы более высокого порядка (как схема Рунге-Кутта), но они потребуют условие сходимости и это условие весьма жесткое для маленьких частиц (потому что та ^ 0).

1.2.2 ГЕ^ЛР-1СЕ

В отличии от LEWICE и ОМЕЯЛ, FENSAP-ICE использует подход Эйлера к описанию движения переохлажденных капель. Такая модель была представлена в 1999г в статье [22]. В ней используются законы сохранения массы и импульса капель:

да

-¡ц+У. (аиа) = 0,

диа С0 Яеа

+ ий.Чий =

д

24К

(иа — иа) + {1—^)7129'

(8)

где:

а - объёмная доля воды,

ыа - безразмерная скорость капель (К/К»),

ыа - безразмерная скорость воздуха,

рwd|иa-иd| „

Ке а =--число Рейнольдса относительно диаметра капель,

¡1

тг рWd^V<Ю

К = ---параметр инерции,

18Ь Л

Са - коэффициент сопротивления капли:

а

V

Са=—(1 + 0.015Re0c[687)lесяиRed > 1300; ^ с

Сс = 0.4, если Red > 1300.

Рг = - число Фруда,

Кю - скорость набегающего потока воздуха, Ь - характерная длина, й - диаметр капли.

Для решения системы водности выбран метод конечного элемента с добавлением стабилизации Петрова-Галеркина.

1.3 Моделирование образования льда

1.3.1 Модель Мессингера

Одной из наиболее распространенных моделей расчета нарастания льда является модель Мессингера 1953 года [23]. Она используется в таких модулях моделирования обледенения , как LEWICE [10], ОМЕЯАЗВ [11], САМОЕ [24]. Эта модель была сформулирована путем применения баланса массы и энергии в конечном контрольном объеме (ячейка, панель) на поверхности (рисунок 1.2).

Рисунок 1.2 - Физические процессы на поверхности [23]

Уравнения баланса массы и энергии записываются в форме:

mimp + m in mice + mevap + mout'

. _ . . P . . (9)

Q imp + Q in Q ice + Q evap + Q out + Q conv + Q rad ■

Изменение массы воды в объеме вычисляется как сумма потока капель mimp и потоков втекающей воды из соседних объемов min, при этом часть ее может превратить в лёд mice, в пар mevap или вытекает из объема mout. В балансе энергии, можно разделить поток тепла на две группы: охлаждение и нагревание. Соответственно получается потоки тепла втекания и вытекания, кроме того есть Qconv - поток тепла конвекция и Qrad - поток тепла излучения. Эти потоки тепла можно записать в следующем виде:

где: f =

Q imp

Q in

Q ice Q evap Q out

Q conv

Q rad

Ш ice

^ imp+m in

m,

imp

mm

ur

cPtW(Td - 27315) + -f

^in[Cp.wiTin -273.15)],

= f (mimp + min) [cp,ice(Tsurf - 273.15) - f],

mevap [^-p,w(Tsurf 273.15) + ^evap],

= [(1 - f)(m^imp + ii in) - m evap ]Cp,w(Tsurf - 273.15),

(10)

Tsurf Te

2

Tcue

2

p,a

0£(Tstrf - О

доля замерзающей воды,

cp,w, cp,iсe, cp,a - соответственно удельная теплоемкость воды, льда, воздуха, Td,Tsurf,Te,Tin - температура капли, на поверхности, на внешней границе граничного слоя, втекающего потока, Levap, Lfus - удельная теплота испарения и плавления льда,

ud, ue, гс - скорость капли, скорость воздуха на внешней границе граничного слоя, фактор восстановления,

о, £ - постоянная Стефана-Больцмана и коэффициент излучения. Для решения этих уравнений, т.е. для определения состояния каждой ячейки (лед, вода, лед-вода), часто используется метод проб и ошибок (например в работе [25]). Процесс решения начинается в точке торможения, где нет втекающих потоков (m.in = 0), которую приходится определять в начале расчета.

1.3.2 Мелководная модель SWIM

В 2000 году была разработана новая модель нарастания льда для программы FENSAP-ICE [26, 27], которая называется «мелководная модель нарастания льда» (Shallow-Water Icing Model или SWIM). Эта модель позволяет предсказывать процесс нарастания льда и течение воды на поверхности.

Рассмотрим плёнку воды на поверхности. На рисунке 1.3 представлены процессы обмена тепла и массы в конечном объеме между твердой поверхностью

и поверхностью пленки воды. Пусть скорость водяной пленки есть функция от координат х = (х1, х2) на поверхности и у- нормали к поверхности:

йт(х,у) = fQt-wa.ll' У),

Рисунок 1.3 - Контрольный объем [26]

где I - напряжения трения воздуха на поверхности, главная движущая сила пленки воды. Поскольку толщина пленки воды маленькая, для простоты предположим, что профиль скорости линейный. Скорость на стенке равна нулю:

Осредняя по толщине пленки, получаем среднюю скорость:

1 Снг кг

йГ (х)=1Т I й}(Х' У) аУ = ^аи (х).

(11)

1 Скг

V

Применяя закон сохранения массы и энергии в контрольным объеме, получаем систему дифференциальных уравнений:

Pw

дк

г

д

+ У(икг)

^еуар,

Pw

<

где:

дкгСр^Т

д

+ У(икгСр^Т)

^^тр + ^^се ^2еуар

сопу чеуар чгай'

(?£тр

Щтр = итьшср,

М2'

+

2

Ш-1тр,

Сеуар 0•5(^evap +

evap^

О.сопу Кс(Т Тгес),

Стай = ^(Т4 - Т^).

Массовой поток в процессе испарения воды и сублимации льда определяется по формуле:

Wlevap

0.7К с у у ¿иг, —

-р,а

где давление насыщенного пара вычисляется по формуле:

ру = 3386(0.0039 + 6.8096 X 10-6Ту2 + 3.5579 X 10-7Ту3), Ту = 72 + 1.8(Т- 273.15). Скорость движения слоя воды вычисляется по формуле:

(14)

(15)

К

и Тщаи

(16)

В системе два уравнения, но три переменные к, Т, т;се, соответственно, толщина, температура пленки воды (в градусах Цельсия) и массовый поток нарастающего льда. Система уравнений является незамкнутой. Для её замыкания добавлены соотношения:

К, > 0, ЩСе > 0, Л,Т > 0, Ш^Т < 0. (17)

Первое соотношение естественно для толщины пленки воды. Второе соотношение

предполагает, что нет плавления ледяной поверхности. Третье и четвертое

означают, что при температуре больше нуля в контрольном объеме есть только

вода, а новый лёд не появляется, в противоположном случае вся вода в объеме

превращается в лёд, а когда температура равна нулю то в объеме одновременно

существуют и вода и лёд.

1.4 Итог обзора литературы

С развитием вычислительной техники, моделирование вязкого течения воздуха становится стандартным. Использование модели вязкого турбулентного газа оказывается более предпочтительным, чем модели невязкого течения с «намазыванием» пограничного слоя.

При моделировании движения переохлажденных капель, по сравнению с подходом Лагранжа, у подхода Эйлера имеется ряд преимуществ: совместимость с солвером движения воздуха, в котором обычно используется подход Эйлера; более простое применение для сложных геометрий; меньшие требования к вычислительным ресурсам, особенно в больших 3D задачах.

Для термодинамического процесса нарастания льда, система SWIM более точно описывает движение пленки воды на поверхности. Для ее решения не надо определить точку торможения, что сложно делать в трехмерных задачах.

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

• уравнения RANS замкнутые моделью турбулентности Спаларта-Аллмараса с учетом шероховатости,

• модель Эйлера движения переохлажденных капель,

• мелководная модель нарастания льда SWIM.

Для решения всех моделей используется метод конечного объёма на неструктурированной сетке.

Глава 2. Математическая постановка и численный метод

2.1 Метод конечного объёма

Рассмотрим закон сохранения некоторой физической величины (масса, импульс, энергия) в контрольном фиксированном объеме V, ограниченном поверхностью

Список литературы диссертационного исследования кандидат наук Нгуен Нгок Шанг, 2022 год

Список литературы

1. FlowVision 3.12.03, User's guide // . - 2021

2. ANSYS 14.5, Apabilities Brochure // ANSYS, Inc. - 2013

3. NUMECA, NUMECA Capability Chart // NUMECA. - 2019

4. COMSOL 5.5, COMSOL Multiphysics Reference Manual // COMSOL. - 2019

5. Autodesk, Autodesk Simulation CFD 2012 // Autodesk. - 2012

6. OpenFoam, User Guide version 9 // OpenFoam. - 2021

7. Thomas D. E., Francisco P., Sean R. C., Trent W. L., Juan J. A., SU2: An Open-Source Suite for Multiphysics Simulation and Design // AIAA Journal. - 2016

8. A. Aksenov, P. M. Byvaltsev, S. V. Zhluktov, K. E. Sorokin, A. A. Babulin, and V. I. Shevyakov, Numerical simulation of ice accretion on airplane surface // AIP Conference Proceedings 2125. - 2019

9. Beaugendre H., Morency F., Habashi W.G, FENSAP-ICE's Three-Dimensional InFlight Ice Accretion Module: ICE3D // Journal of Aircraft, Vol. 40, No. 2.. - 2003

10. Ruff G., Berkowitz B., Users Manual for the NASA Lewis Ice Accretion Prediction Code (LEWICE) // NASA Technicall Report. - 1990

11. Hedde T., Guffond D., ONERA Three-Dimensional Icing Model // AIAA Journal, Vol. 33, No. 6. - 1995

12. Bourgault-Cote S., Laurendeau E. , Two-Dimensional/Infinite Swept Wing Ice Accretion Model // AIAA paper 2015-0535. - 2015

13. Босняков С. М., Власенко В. В., Енгулатова М. Ф., Зленко Н. А., Матяш С. В., Михайлов С. В., Программный комплекс EWN // Свидетельство о государственной регистрации программы для ЭВМ № 2008610227. - 2008

14. Wright W.B., Gent R.W., GuffondD, DRA/NASA/ONERA Collaboration on Icing Research. Part II-Prediction of Airfoil Ice Accretion // NASA Technicall Report. -1997

15. Trontin P., Kontogianis A., Blanchard G., Villedieu P., Description and assessment of the new ONERA 2D icing suite IGLOO2D // AIAA 2017-3417. - 2017

16. Beaugendre H., Morency F., Habashi W.G., Development of a Second Generation In-Flight Icing Simulation Code // Journal of Fluids Engineering, Vol. 128 / 387. -2006

17. Geuzaine C., Remacle J., Gmsh Reference Manual // Gmsh 4.9.3. - 2022

18. SALOME, User's Documentation // SALOME platform. - 2021

19. Colin S.B, Mark G. P., Users Manual for the NASA Lewis Three-Dimensional Ice Accretion Code (LEWICE 3D) // NASA Technical Memorandum 105974. - 1993

20. Aupoix B., Spalart P. R., Extensions of the Spalart-Allmaras turbulence model to account for wall roughness // International Journal of Heat and Fluid Flows, Vol. 24. -2003

21. William B. W., User manual for the NASA Glenn Ice Accretion Code LEWICE Version 2.2.2 // NASA/CR -- 2002-211793. - 2002

22. Bourgault Y., Habashi W. G, Dompierre J., Baruzzi G. S., A Finite Element Method Study of Eulerian Droplets Impingement Models // International Journal for Numerical Methods in Fluids, Vol. 29. - 1999

23. Messinger B.L., Equilibrium Temperature of an Unheated Icing Sureface as a Function of Air Speed // Journal of the Aeronautical Sciences, Vol. 20, pp. 29-42. -1953

24. Morency F., Tezok F., Anti-Icing System Simulation Using CANICE // Journal of Aircraft. - 1999

25. Зыонг Де Тай, Обтекание планера гражданского самолета в условиях начальной стадии обледенения // Московский Физико-Технический Институт. -2018

26. Bourgault Y., Beaugendre H., Habashi W. G., Development ofa Shallow Water Icing Model in FENSAP-ICE // Journal of Aircraft, Vol. 37. - 2000

27. Beaugendre H., A PDE-Based 3D Approach to In-Flight Ice Accretion // McGill Univeristy. - 2003

28. Shima E., Kitamura K., Haga T., Green-Gauss/Weighted-Least-Squares Hybrid Gradient Reconstruction for Arbitrary Polyhedra Unstructured Grids // AIAA Journal. T.51, № 11.. - 2013

29. Barth T., Jesperson D., The design and application of upwind schemes on unstructured meshes // AIAA Paper. - 1989

30. Michalak K., Ollivier-Gooch C., Limiters for Unstructured Higher-Order Accurate Solutions of the Euler Equations // AIAA Paper. - 2008

31. Venkatakrishnan VOn the accuracy of limiters and convergence to steady state solution // AIAA Paper. - 1993

32. Власенко В., О математическом подходе и принципах построения численных методологий для Пакета Прикладных Программ EWT-ЦАГИ // Труды ЦАГИ, № 2671. - 2007

33. Nikuradse J.., Laws of Flow in Rough Pipes // NASA /TM-1292. - 1937

34. Grigson C., Drag losses of new ships caused by hull finish // J. Ship Res. 36. -1992

35. Shin J., Berkowitz B., Chen H., Cebeci T., Prediction of ice shapes and their effect on airfoil performance // AIAA paper 91-0264. - 1991

36. Aupoix B., Improved heat transfer predictions on rough surfaces // International Journal of Heat and Fluid Flow. - 2015

37. Chedevergne F., Aupoix B., Accounting for wall roughness effects in turbulence models: a wall function approach // EUCASS2017-372. - 2017

38. Власенко В.В., Расчетно-теоретические модели высокоскоростных течений газа с горением и детонацией в каналах // ФГУП ЦАГИ. - 2017

39. S. K. Jung, R. S. Myong, T. H. Cho, Development of Eulerian Droplets Impingement Model Using HLLC Riemann Solver and POD-Based Reduced Order Method // 41st AIAA Fluid Dynamics Conference and Exhibit. - 2011

40. S. Keita, Y. Bourgault, Eulerian models with particle pressure for air-particle flows // European Journal of Mechanics / B Fluids. - 2019

41. E. F. Toro, Riemann Solvers and Numerical Methods for Fluid Dynamics // Springer-Verlag Berlin Heidelberg. - 2009

42. Langmuir I., Blodgett K. B., A mathematical investigation of water droplet trajectories // Army Air Forces Technical Report. - 1946

43. Iuliano E., Montreuil E., Norde E., Van der Weide E.T.A, Hoeijmakers H.W.M, Modelling of Non-Spherical Particle Evolution for Ice Crystals Simulation with an Eulerian Approach // SAE Technical Paper. - 2015

44. Ellen Norde, Eulerian Method for Ice Crystal Icing in Turbofan Engines // University of Twente. - 2017

45. Ganser G. H., A rational approach to drag prediction of spherical and nonspherical particles // Power Technology 77(2):143-152. - 1993

46. Lavoie P., Bourgault-Cote S., Laurendeau E., Numerical algorithms for ifinite swept wing ice accretion // J. Computers and Fluids. - 2018

47. S. Keita, Eulerian Droplet Model: Mathematical analysis, improvement and applications // University of Ottawa. - 2018

48. Sana Keita, Yves Bourgault, Eulerian droplet model: Delta-shock waves and solution of the Riemann problem // Journal of Mathematical Analysis and Applications. - 2019

49. F Bouchut, On zero pressure gas dynamics // Advances in Kinetic Theory and Computing. - 1994

50. Papadakis M., Rachman A., Wong S.C., Yeong H.W., Hung K.E., Vu G.T., Bidwell C.S., Water droplet impingement on simulated glaze, mixed, and rime ice accretions // NASA/TM-2007-213961. - 2007

51. Hosni M., Coleman H., Taylor R., Measurements and calculations of rough-wall heat transfer in the turbulent boundary layer // J. Heat Mass Transf. Vol. 34. - 1991

52. Hosni M., Coleman H., Garner J., Taylor R., Roughness element shape effects on heat transfer and skin friction in rough-wall turbulent boundary layers // J. Heat Mass Transf. Vol. 36. - 1993

53. Olumide F. O., Upgrading and Qualification of a Turbulent Heat Transfer Test Facility // Mississippi State University, Thesis. - 2002

54. Healzer J., Moffat R., Kays W., The turbulent boundary layer on a porous, rough plate: Experimental heat transfer with uniform blowing // AIAA/ASME 1974 Thermophysics and Heat Transfer Conference. - 1974

55. Airbus Industrie, Getting to grips with cold weather operations: A Flight Operation View // AIRBUS INDUSTRIE. - 2000

56. Yannick L. B., Jean-Marc G., Magneto-Elastic Resonance: Principles, Modeling and Applications // IntechOpen. - 2017

57. Шевяков В. И., Разработка теоритических основ и практических методов реализации аэродинамического совершенства самолётов транспортной категории с учётом выполнения сертификационных требований по безопасности полёта // МГТУ ГА, Диссертация. - 2017

58. ParaView 5.6, The ParaView Guide // Kitware Inc. - 2019

59. Шевяков В.И., К вопросу обеспечения безопасности полетов в условиях обледенения // Научный вестник МГТУ ГА, 2011, No 172. - 2011

60. Бобцов В.А., Зленко Н.А., Автоматизированная диалоговая система аппроксимации экспериментальных данных // Труды ЦАГИ, No2522. - 1993

61. Зленко Н.А., Программа анализа и аппроксимации дискретно заданных одномерныхзависимостей (APPEX) // Свидетельство об официальной регистрации программыдля ЭВМ No2007614013. - 2007

62. Дрейпер Н., Смит Г., Прикладной регрессионный анализ // Москва "Финансы и Статистика. - 1973

63. Адлер Ю.П., Маркова Е.В., Грановский Ю.В, Планирование эксперимента при поискеоптимальных условий // Москва "Наука". - 1976

Отпечатано с оригинал-макетов Заказчика в типографии "Переплетофф" Адрес: г. Долгопрудный, ул. Циолковского, 4. Тел: [phone]. www.perepletoff.ru Формат 210 х 297 мм. Бумага офсетная. Печать цифровая. Тираж 11 экз. Твердый переплет. Заказ № . 28.03.22 г.