УДК 533
БОТ: 10.15350/17270529.2019.3.43
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ И ВЫЧИСЛИТЕЛЬНЫЙ ЭКСПЕРИМЕНТ ДЛЯ ОПРЕДЕЛЕНИЯ ПАРАМЕТРОВ ТЕПЛОМАССОПЕРЕНОСА ТЕРМОВИБРАЦИОННОГО КОНВЕКТИВНОГО РЕЖИМА ТЕЧЕНИЯ В ТОПЛИВНОМ БАКЕ ЖРД
1МОРМУЛЬ Р. В., 2сальников а. ф., 3павлов д. а., 1гладков а. н.
1 Пермский военный институт войск национальной гвардии РФ, 614112, г. Пермь, ул. Гремячий Лог, 1
2 Пермский национальный исследовательский политехнический университет, 614990, г. Пермь, Комсомольский пр., 29
ПАО НПО Искра, 614038, г. Пермь, ул. Академика Веденеева, д. 28
АННОТАЦИЯ. Численно исследована структура осреднённого и пульсационного течения в топливном баке жидкостного ракетного двигателя (ЖРД) при экзотермических и обратных эндотермических процессах горения. Изучены механизмы трансформации структуры течения в зависимости от ориентации векторов напряженности гравитационного и вибрационного полей. Проведен анализ степени влияния параметров управления и контроля конвективными процессами неизотермических продуктов сгорания ЖРД на распределение характеристик тепломассопереноса в переменных по величине и направлению силовых полях в условиях пониженной гравитации. Получены численные оценки по распределению кинетической энергии осредненного и пульсационного движения конвективной системы.
КЛЮЧЕВЫЕ СЛОВА: термовибрационная конвекция, ЖРД, численное моделирование, тепломассоперенос, метод конечных разностей, полиномы Чебышева, дискретное преобразование Лапласа.
Устойчивый интерес исследований, посвященный изучению конвективного теплообмена в элементах двигательных установок на жидком топливе при длительном космическом полете под воздействием разного рода вибраций, обусловлен развитием космических технологий.
Топливные баки, используемые в жидкостных ракетных двигательных установках (ДУ), в общих случаях могут быть предназначены для:
- хранения основного запаса жидкого ракетного топлива, т.е. топлива, предназначенного для получения всей или основной доли тяги двигателя. Эти емкости называют основными топливными баками ДУ;
- хранения вспомогательного жидкого ракетного топлива, предназначенного для работы вспомогательных систем и агрегатов ДУ;
- выполнения роли промежуточных (буферных) емкостей, предназначенных для периодического накопления и расходования топлива в процессе работы двигателя.
Конструкция баков должна быть устойчива по отношению к воздействиям внешних (аэродинамических) и внутренних (давление вытеснения топлива, статическое и динами -ческое воздействия массы топлива, находящегося в баках) сил, инерционных перегрузок и должна обеспечить максимальную полноту выработки топлива, а также минимальное растворение (или конденсацию) газа, находящегося в баке в компоненте топлива.
При хранении топлива наиболее важными характеристиками являются рост давления в баке и температуры на входе в двигатель, допустимой из условия появления кавитации. Поэтому представляет интерес не только теплопередача, но и температура на межфазной поверхности раздела, для расчета которой необходимо привлекать более полную модель перемешивания жидкости в баке с учетом положения межфазной границы в объеме бака. Определение последней в условиях невесомости представляет самостоятельную задачу. В предельном случае при полете частично заполненного топливом бака в режиме теоретической невесомости (§ ~ 0,1) межфазная поверхность, в случае топлива, смачивающего стенки, сворачивается в пузырь, который может в зависимости от возмущений располагаться в различных частях бака.
Одним из возможных и наиболее простых расчетных случаев является его расположение внутри, в центре бака. При этом единственным механизмом перемешивания в режиме теоретической невесомости является термокапиллярная конвекция Марангони, которая, как и тепловая гравитационная конвекция, проявляется при потере устойчивости механического равновесия или при его отсутствии.
В данной статье в высокочастотном приближении сформулированы уравнения термовибрационной конвекции в топливном баке ЖРД, описывающие медленное осредненное конвективное движение, возникающее на фоне быстрых пульсаций полей скорости, давления и температуры.
При описании термовибрационной конвекции, как правило, используются уравнения Зеньковской-Симоненко [1]. Систематическое описание основных положений термовибрационной конвекции представлено в монографии [2]. Процесс внутреннего конвективного движения газа сопровождается тепловыделением, обусловленным протеканием в реагирующей среде химических реакций [3].
На рис. 1 представлена геометрия топливного бака ЖРД, введены обозначения:
(х, у, z) - декартовая система координат; п, к - единичный вектор нормали и вибраций; составляющие углы а, /, соответственно с осью х; д - напряженность гравитационного поля.
Для рассмотрения поведения газообразных продуктов течения используются уравнения тепловой вибрационной конвекции [4], хорошо обоснованные, как теоретически [5 - 8], так и экспериментально [9, 10].
Согласно принципу суперпозиции, скорость, температура и давление в системе координат, связанной с полостью, представляются в виде сумм двух слагаемых V + г), Т и р + 5, где V, Т и р - осредненные величины, медленно меняющиеся во времени, а т),^ и 5- пульсационные величины, осциллирующие с частотой с .
Рис. 1. Геометрия задачи термовибрационной конвекции
Конвективное движение продуктов сгорания ЖРД сопровождается энерговыделением, обусловленным протеканием экзо- и эндотермической реакции в форме Аррениуса - Франк-Каменецкого:
[вг = во
1 ^ (1) б2 = вое т (1 -Л)
В соотношении (1) введены обозначения: б - начальная мощность внутренних источников тепла; в , в - теплоты образования при экзотермических и обратных эндотермических химических реакциях, соответственно; - параметр энерговыделения; Л - параметр поглощения энерговыделения при обратных эндотермических процессах.
Систему уравнений термовибрационной конвекции продуктов сгорания ЖРД запишем в безразмерной форме следующим образом:
+ = Аф + игЧТ х п + (VW • k) хЧТ
дф
дТ 1 Fr
— + V•ЧТ = — АТ + —О ,I = 1,2 дt Рг Рг г
Ащ + ф =0 АР + ЧТ х к = 0
V = Чхщ W = Чх Р
V •V = 0
V ^ = 0
(2)
В системе уравнений (2) введены: 1) безразмерные критерии подобия:
Р г , г
иг =----число Грасгофа; иг =
У X 2
1 ( аЬаОоР
УХ
— вибрационное число Грасгофа;
кЕх ехр
'- Ех ^
Рк =
ЕР О
0
У
Е
— число Франк-Каменецкого; Рг = — — число Прандтля.
X
2) геометрические параметры модели:
Е — радиус внутреннего цилиндра; Е2 = Е + Р — радиус внешнего цилиндра; Р — величина зазора.
3) соленоидальные вектор-функции:
V — поле осредненной скорости конвективного течения;
W — амплитуда скорости пульсационного течения; Щ,ф — функция тока и вихря скорости осредненного течения;
Р — функция тока пульсационного движения газообразных продуктов сгорания.
4) дифференциальные операторы цилиндрической системы координат:
^ г д А д2 1 д 1 д2 д2
V = е —т — оператор Гамильтона; А = —у +----+ -у —у +--у — оператор Лапласа.
д^ дг г дг г д$ дг
5) единичные вектора нормали п и вибраций к.
6) Е — энергия активации химической реакции; Е — универсальная газовая постоянная, Т — температура газа.
Из условия прилипания осредненной скорости течения, условия непротекания пульсационной составляющей формируются граничные условия для функций тока:
= 0. (3)
Границы полости предполагаются изотермическими:
дщ н Щ = —!— = Р\
дг
Т
Г=Е
= 800 К;
Т
г=К.
= 50 К;
Т\,=0 = Т,
(4)
<
МЕТОД РЕШЕНИЯ
В качестве численного решения конвективной граничной задачи (1) - (4) использован метод конечных разностей [11].
Для записи уравнений (2) в конечно-разностной форме введем обозначения: г = IАг, Зк = кА$\ где , = 0,1,...,I; к = 0,1,..., К.
Аг = 1/1, АЗ* = 2п/ К, АЗ = 28т(АЗ*/2). (5)
Для функций и конечно-разностных аппроксимаций дифференциальных операторов в узлах сетки применим общепринятые обозначения:
/к = /(г,З), / = (/+1 -/Л/2аГ,
/з = (Л+1 - Л-1)/2АЗ, (6)
/з = (Л+1 - 2/ + Лк-1)/2АЗ
Лгг = (/+1 - 2Л + Л^1)/2Аг,
Лз = (/+1,к+1 - Л-1,к+1 - Л+1,к-1 + Л-1 , к-1)/ 4АгА3,
/ = С/>+1 - /-1)/2Аг
Для определения равновесной температуры в газе использован операторный метод на базе дискретного преобразования Лапласа [12].
Операторная форма записи уравнения теплопроводности с учетом логарифмирования примет следующий вид:
а) для случая экзотермической реакции
1п (уТ (г, 5) - Т0 ) = 1п
( {д 2Т (г, у) 1 дТ (г, у) ^ Ек
\
Рг-1
V V
дг2
г дг
+ — 0>еАЕТ ) Рг
(7)
у
б) для случая эндотермической реакции
1п (уТ (г, у) - Т ) = 1п
Рг-1
гд Т (г, у) 1 дТ (г, у)^
дг2
г
дг
Ек
у
+ — б0еТ )(1 -Л)
Рг
(8)
На рис. 2 представлены профили равновесной температуры в газе для различных значений числа Франк-Каменецкого и параметров энерговыделения конвективной системы.
Я, м
1.0
1.5
2.0
Рис. 2. Профили равновесной температуры в газе (Рг и 0,78) для различных значений числа Франк-Каменецкого и параметров энерговыделения
1: Ек = 0,6; 2: Ек = Л = 0,4; 3: Ек = 0,3
Максимум равновесной температуры в газе достигается в случае 1 рис. 2.
А
Е
1) Свободный конвективный режим
Конвективное движение неравновесной среды происходит под действием статического поля силы тяжести, и обусловлено температурной неоднородностью поля температуры, определяемого условиями подогрева на границах топливного бака.
На рис. 3, а визуализированы линии тока конвективной системы при протекании экзотермических реакций, на рис. 3, б — при обратных эндотермических химических реакциях.
В полости формируется двухъярусная четырехвихревая серповидная конвективная структура. Масштаб вихря определяется геометрическими размерами цилиндрического канала топливного бака.
а) б)
Рис. 3. Карта функции тока при свободном конвективном режиме:
а) - иг = 1000; Я, = 0,2 м; Щ = 2,2 м; Fk = 0,4; Grv = р = 0; а = ж/2
б) - иг = 1500; Е = 1,2м; Я2 = 2,2м; 2 = Бк = 0,6; Grv =р = 0; а = ж/2
2) Термовибрационный конвективный режим
При запуске и работе ЖРД требуется обеспечить гарантированный забор жидких компонентов из топливных баков как в условиях невесомости, так и при наличии переменных по величине и направлению ускорений.
Доминирующее влияние вибраций на конвекцию в условиях пониженной гравитации также определяет практическую необходимость исследований, результаты которых необходимо использовать в области космических технологий, активно развиваемых в настоящее время. Отметим, что высокая стоимость экспериментов по определению параметров тепломассопереноса в топливных баках ЖРД, проводимых в башнях невесомости, даёт почву для поиска альтернативных методик исследования, одной из которых является численное моделирование.
При доминировании вибрационного поля над напряженностью гравитационного, вихревая структура нижнего яруса вытесняет структуру верхнего (рис. 4). В случае взаимно ортогональной ориентации силовых полей (гравитационного и вибрационного) осевая симметрия течения сохраняется, в общем случае (произвольной ориентации) — нарушается.
Механизм трансформации вихревых структур течений термовибрационного режима объясняется преобладанием кинетической энергией осреднённого движения вибрационного конвективного режима течения над кинетической энергией осреднённого течения при реализации свободного.
Рис. 4. Карта функции тока осредненного течения при термовибрационном конвективном режиме:
а) - иг = 100; Е = 1,2 м; Я2 = 2,2 м; Бк = 0,4; Grv = 1000; / = 1,75ж; а = ж/6; ? = 5 с
б) - иг = 10; Е = 0,2 м; Щ = 2,8 м; 2 = Бк = 0,6; Grv = 1500; / = ж/10; а = ж/3; t = 5 с
в) - иг = 50; Е = 0,2м; Щ = 2,2м; 2 = Бк = 0,6; Grv = 1500; / = ж; а = ж/2; t = 5с
г) - иг = 50; Е = 1,2м; Я2 = 2,2м; Бк = 0,4; Grv = 1500; / = ж; а = ж/2; t = 5 с
Таким образом, при увеличении происходит относительное уменьшение
интенсивности внешних вихрей и увеличение интенсивности внутренних, увеличивается вклад кинетической энергии пульсационного движения по сравнению с кинетической энергией среднего течения:
2п 2п
-^Ц | Пс| 2гс1гМ>-^1 | I гс I 2гсгг£^т9 . (9)
0 Й! О
Данный процесс сопровождается перераспределением тепловых потоков на поверхности коаксиальных цилиндров топливного бака и смещением центров вихрей конвективной системы.
Пульсационные характеристики среды определяются по осредненным согласно формулам:
f = ■ ЧТ) sin ílt
rj = 7 2 Grv Же o sil t (10)
ii = = S t - безразмерная частота колебаний (параметр Стокса).
Численное моделирование термодинамических процессов, происходящих в реальных конструкциях, в условиях невесомости и неопределённости граничного и фазового состояния жидкого топлива представляет собой сложнейшую задачу.
Наиболее точные результаты получаются при проведении физического моделирования на конструктивно подобных моделях, в которых соблюдается подобие не только геометрических, но и динамических и других параметров процесса.
Для конструктивно подобных моделей необходимо обеспечить не только геометрические параметры, но и габаритно-массовые, жесткостные и другие характеристики топливной системы.
Карты распределения функции тока пульсационного течения представлены на рис. 5.
Рис. 5. Карты функции тока пульсационного течения
а) - Gr = 100; Grv = 1000; R = 0,5м; R = 1,5м; Fr = 0,6; а = п/6; / = 2^/3; t = 5с
б) - Gr = 100; Grv = 2000; R = 0,2м; R = 1,2м; Л = Fr = 0,4; а = п/6; /3 = 5ж/9; t = 7с
ВЕРИФИКАЦИЯ
В целях верификации, решение (1) - (5), полученное методом конечных разностей сравнивалось с численным решением, реализованным псевдо-спектральным методом с использованием полиномов Чебышева и разложением Фурье [8].
Введем в рассмотрение новую радиальную координату /л = 2(г - 0,5(А +1)) при условии, что г е [A /2, A /2 +1]. Таким образом, спектральное представление функции тока и поля температуры осредненного течения примет вид:
M/2 N ~
Г(М,3, о = £ Ё атп (ОТ ^)е""3 (11)
m=-M /2+1 п=0 М/2 N
Т(р,3, ^ = Ё Ё Ьтп (t )Т~п (0е'т3 (12)
т=-М /2+1 п=0
В соотношениях (11) - (12): Тп(t) = (1 -ц2)ТпV); ~(t) = (1 - V2)2Тп(ц).
Результаты численного моделирования термовибрационного конвективного режима представлены на рис. 6.
в) г)
Рис. 6. Карта функции тока осредненного течения при термовибрационном конвективном режиме:
иг = 100; Е = 0,05 м; Щ = 1,05м; Бк = 0,4; Grv = 2500;/ = 2ж/3; а = ж/6; t = 8с а) — Метод конечных разностей
б) — псевдо-спектральный метод, М = N = 18
в) — псевдо-спектральный метод, М = N = 28
г) — псевдо-спектральный метод, М = N = 100
Следует отметить, что решение термовибрационной конвективной задачи, полученное псевдо-спектральным методом хорошо согласуется с численным решением, реализованным методом конечных разностей (рис. 6, а и 6, г). Несмотря на то, что гармоники в разложении Фурье являются собственными функциями оператора Лапласа в уравнениях Пуассона, скорость сходимости метода конечных разностей на 26 % выше, чем у псевдо-спектрального.
ЗАКЛЮЧЕНИЕ
Разработанная математическая модель и комплекс вычислительных алгоритмов позволяет адекватно исследовать структуру конвективных течений в условиях сильной гравитации, а также в переменных по величине и направлению вибрационных полей в условиях микрогравитации при модуляциях гравитационного поля (инерционные силы и перегрузки). Данный подход численного моделирования позволяет осуществлять контроль и управление конвективными процессами неизотермических продуктов сгорания на распределение характеристик тепломассопереноса в топливных баках ЖРД.
СПИСОК ЛИТЕРАТУРЫ
1. Зеньковская С. М. Исследование конвекции в слое жидкости при наличии вибрационных сил // Известия АН СССР. Механика жидкости и газа. 1968. № 1. С. 55-58.
2. Gershuni G. Z., Lyubimov D. V. Thermal Vibrational Convection. John Wiley & Sons, 1998. 358 p.
3. Франк-Каменецкий Д. А. Диффузия и теплопередача в химической кинетике. М.: Наука, 1987. 502 с.
4. Повицкий А. С., Любин Л. Я. Основы динамики и тепломассообмена жидкостей и газов при невесомости. М.: Машиностроение, 1972. 252 с.
5. Полежаев В. И., Соболева Е. Б. Тепловая гравитационная и вибрационная конвекция околокритического газа в условиях микрогравитации // Известия Российской академии наук. Механика жидкости и газа. 2000. № 3. С. 70-80.
6. Бабский В. Г., Копачевский Н. Д., Мышкис Л. Д. и др. Гидродинамика невесомости. М: Наука, 1976. 504 с.
7. Бурдэ Г. И. Численное исследование конвекции, возникающей в модулированном поле внешних сил // Известия АН СССР. Механика жидкости и газа. 1970. Т. 5, № 2. С. 196-201.
8. Острах С. Роль конвекции в технологических процессах, проводимых в условиях микрогравитации // В сб. статей: Космическая технология / под ред. Л. Стега; пер. с англ. под. ред. А. С. Охотина. М.: Мир, 1980. С. 9-27.
9. Зюзгин А. В., Иванов А. И., Полежаев В. И., Путин Г. Ф., Соболева Е.Б. Исследование околокритической жидкости в условиях микрогравитации: эксперименты на станции «Мир» и численное моделирование // Космонавтика и ракетостроение. 2000. № 19. C. 56-63.
10. Иванова А. А., Козлов В. Г. Экспериментальное исследование тепломассопереноса в условиях вибрационной конвекции // Тезисы докладов III Всесоюзного семинара по гидромеханике и тепломассопереносу в невесомости. Черноголовка, ИФТТ АН СССР, 1984.
11. Ладыженская О. А. Метод конечных разностей в теории уравнений с частными производными // Успехи математических наук. 1957. Т. 12, № 5(77). С. 123-148.
12. Сергиенко А. Б. Цифровая обработка сигналов. Учебник для вузов. 2-е изд. СПб.: Питер, 2006. 751 с.
MATHEMATICAL SIMULATION AND NUMERICAL EXPERIMENT FOR DETERMINING THE PARAMETERS OF HEAT AND MASS TRANSFER OF THERMAL VIBRATIONAL CONVECTIVE FLOW IN THE FUEL TANK OF LRM
'Mormul R. V., 2Salnikov A. F., 3Pavlov D. A., 'Gladkov A. N.
1 Perm Military Institute of National Guard Forces, Perm, Russia
2 Perm National Research Polytechnic University, Perm, Russia
3 PJSC Research and Production Association "Iskra", Perm, Russia
SUMMARY. A numerical study of the structure as averaging and pulsatile flow in the fuel tank of liquid rocket motor (LRM) for exothermic and endothermic reverse combustion processes. The mechanisms of transformation of the structure of flow depending on the orientation of the vectors of gravity and vibration fields. The analysis of the degree of influence of the control parameters and control of non-isothermal processes convective LRM combustion products in the distribution of heat and mass characteristics in varying magnitude and direction of the force fields in low gravity. Numerical estimates are obtained for the distribution of kinetic energy of the averaged and pulsating motion of the convective system.
KEYWORDS: thermal vibrational convection, LRM, numerical simulation, heat and mass transfer, method of finite difference, Chebyshev polynomials, discrete Laplace transform.
REFERENCES
1. Zen'kovskaya S. M. Study of convection in a liquid layer with vibrating forces. Fluid Dynamics, 1968, vol. 3, iss. 1, pp. 35-36. https://doi.org/10.1007/BF01016232
2. Gershuni G. Z., Lyubimov D. V. Thermal Vibrational Convection. John Wiley & Sons, 1998. 358 p.
3. Frank-Kamenetskiy D. A. Diffuziya i teploperedacha v khimicheskoy kinetike [Diffusion and heat transfer in chemical kinetics]. Moscow: Nauka Publ., 1987. 502 p.
4. Povitskiy A. S., Lyubin L. Ya. Osnovy dinamiki i teplomassoobmena zhidkostey i gazov pri nevesomosti [Fundamentals of dynamics and heat and mass transfer of liquids and gases in zero gravity]. Moscow: Mashinostroenie Publ., 1972. 252 p.
5. Polezhaev V. I., Soboleva E. B. Thermo-gravitational and vibrational convection in a near-critical gas in microgravity. Fluid Dynamics, 2000, vol. 35, iss. 3, pp. 371-379. https://doi.org/10.1007/BF02697750
6. Babskiy V. G., Kopachevskiy N. D., Myshkis L. D. i dr. Gidrodinamika nevesomosti [Hydrodynamics of weightlessness]. Moscow: Nauka Publ., 1976. 504 p.
7. Burde G. I. Numerical study of the onset of convection in a modulated external force field. Fluid Dynamics, 1970, vol. 5, iss. 2, pp. 344-349. https://doi.org/10.1007/BF01080256
8. Ostrakh S. Rol' konvektsii tekhnologicheskikh protsessakh, provodimykh v usloviyakh mikrogravitatsii [The role of convection in technological processes carried out under microgravity]. V sb. statey: Kosmicheskaya tekhnologiya [In coll. Articles: Space Technology]. Pod red. L. Stega; per. s angl. pod. red. A. S. Okhotina. Moscow: Mir Publ., 1980, pp. 9-27.
9. Zyuzgin A. V., Ivanov A. I., Polezhaev V. I., Putin G. F., Soboleva E. B. Issledovanie okolokriticheskoy zhidkosti v usloviyakh mikrogravitatsii: eksperimenty na stantsii «Mir» i chislennoe modelirovanie [Investigation of near-critical fluids in microgravity: experiments at the station "Mir" and numerical simulation]. Kosmonavtika i raketostroenie [Space and rocket science], 2000, no. 19. C. 56-63.
10. Ivanova A. A., Kozlov V. G. Eksperimental'noe issledovanie teplomassoperenosa v usloviyakh vibratsionnoy konvektsii [Experimental study of heat and mass transfer in a vibration convection]. Tezisy dokladov III Vsesoyuznogo seminara po gidromekhanike i teplomassoperenosu v nevesomosti [Abstracts of the III All-Union seminar on hydromechanics and heat and mass transfer in zero gravity]. Chernogolovka, IFTT AN SSSR, 1984.
11. Ladyzhenskaya O.A. Metod konechnykh raznostey v teorii uravneniy c chastnymi proizvodnymi [The method of finite differences in the theory of partial differential equations]. Uspekhi matematicheskikh nauk [Advances in Mathematical Sciences], 1957, vol. 12, no. 5(77), pp. 123-148.
12. Sergienko A. B. Tsifrovaya obrabotka signalov [Digital signal processing]. Uchebnik dlya vuzov. 2-e izd. St. Petersburg: Piter Publ., 2006, 751 p.
Мормуль Роман Викторович, доцент кафедры «ВМКСиС» ПВИ ВНГ РФ, e-mail: rmormul@yandex. ru
Сальников Алексей Федорович, доктор технических наук, профессор кафедры «РКТ и ЭС» ПНИПУ, e-mail: afsalnikov1@mail. ru
Павлов Дмитрий Александрович, инженер-конструктор 2 категории ПАО НПО «Искра», e-mail: pal0 7@yandex. ru
Гладков Алексей Николаевич, кандидат технических наук, начальник кафедры «ВМКСиС» ПВИ ВНГ РФ, e-mail: alex-gladkov@mail. ru