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

• СПИСОК СОКРАЩЕНИЙ, УСЛОВНЫХ ОБОЗНАЧЕНИЙ И СИМВОЛОВ 3

• ВВЕДЕНИЕ

  4

• 1. ОБЩИЕ СВЕДЕНИЯ О РАКЕТЕ 3М-14 5
• 2. РАСЧЕТ ДИНАМИКИ ВЫХОДА РАКЕТЫ ИЗ ШАХТНОЙ ПУСКОВОЙ УСТАНОВКИ (аналитический расчет) 8
• 2.1 Исходные данные 8
• 2.2 Постановка задачи 8
• 2.3 Проектирование ШПУ 9
• 2.4 Расчет параметров выхода ракеты при заданной схеме 22
• Выводы 27
• Анализ полученных результатов и заключение 27
• 3. МЕТОДИКА ЧИСЛЕННОГО ЭКСПЕРИМЕНТА 28
• 3.1 Постановка задачи 28
• 3.2 Общие сведения о проведении численного эксперимента газовой динамики 28
• 3.3 Компьютерные пакеты для численного решения задач газовой динамики 37
• 3.4 Общие сведения о методе контрольного объема 39
• Выводы 42
• 4. РАСЧЕТ ДИНАМИКИ ВЫХОДА РАЕТЫ ИЗ ШАХТНОЙ ПУСКОВОЙ УСТАНОВКИ (численный эксперимент) 44
• 4.1 Постановка задачи 44
• 4.2 Стратегия постановки численного эксперимета 44
• 4.3 Анализ полученных результатов 47
• Выводы 49
• 5. ОПРЕДЕЛЕНИЕ АЭРОДИНАМИЧЕСКИХ НАГРУЗОК 50

• Введение

  50

• 5.1 Постановка задачи 52
• 5.2 Подготовка и проведение численного эксперимента 52
• 5.3 Анализ полученных результатов 56
• 5.3.1 Вертикальная составляющая потока () 57
• 5.3.2 Поперечная составляющая потока () 60
• 5.3.3 Сложение составляющих потока (;) 66
• 5.4 Определение моментов 71
• Выводы 79
• ВЫВОДЫ 81

• ЗАКЛЮЧЕНИЕ

  82

• БИБЛИОГРАФИЧЕСКИЙ СПИСОК 83
• ПРИЛОЖЕНИЕ 85

СПИСОК СОКРАЩЕНИЙ, УСЛОВНЫХ ОБОЗНАЧЕНИЙ И СИМВОЛОВ

ШПУ — шахтная пусковая установка;

АПЛ — атомная подводная лодка;

СУ — система управления;

П.П.П. — программный прикладной пакет;

ЭВМ — электронная вычислительная машина;

ПЭВМ — персональная электронная вычислительная машина;

МКО — метод контрольного объема;

ЛА — летательный аппарат;

АДХ — аэродинамические характеристики;

АДК — аэродинамические коэффициенты;

p — давление в подракетном пространстве;

Re — число Рейнольдса;

T_к — температура в камере сгорания;

T₂ — температура в «подракетном» пространстве;

F_x, F_y, F_z — проекция силы на ось X, Y, Z, соответственно;

X — сила лобового сопротивления;

Y — подъемная сила;

Z — боковая сила;

С_x — коэффициент подъемной силы;

С_y — коэффициент нормальной силы;

х_ц.д. — положение центра давления;

С_p — коэффициент давления.

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

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

Перед исполнителем поставлена задача создания схемы расчета динамики выхода ракеты из шахты и определения действующих на нее нагрузок. Решение поставленной задачи ведется на примере ракеты 3М-14.

В найденных литературных источниках тема исследования динамики выхода ракеты из ШПУ недостаточно освещена, ввиду специфичности тематики.

На основе вышеизложенного материала и по мере выполнения задания, в работе сформулированы следующие задачи:

1) Решение задачи динамики выхода ракеты из пусковой установки (газодинамическая задача);

2) Определение А.Д.Х. при поперечном обтекании корпуса ракеты под действием ветровой нагрузки.

1. ОБЩИЕ СВЕДЕНИЯ О РАКЕТЕ 3М-14

Ракета 3М-14 входит в состав системы «Клаб». «Клаб» (Club) — ракетная система с крылатыми ракетами надводного, подводного, наземного и авиационного базирования.

· История создания

Ракетная система «Клаб» (Club) разработана и производится (основные элементы) ОКБ «Новатор» (г. Екатеринбург). Первый испытательный пуск противокорабельной ракеты (ПКР) системы, по данным СМИ, состоялся с атомной подводной лодки (АПЛ) на Северном флоте в марте 2000 г., второй — в июне того же года с дизельной подводной лодки (ДПЛ) проекта 877 Балтийского флота. Оба пуска были признаны успешными.

Первым основным элементом системы является универсальная ракета «Альфа», которая была продемонстрирована в 1993 г. (через 10 лет после начала ее разработки) на выставке вооружений в Абу-Даби и на международном авиакосмическом салоне МАКС-93 в г. Жуковский. В том же году она была принята на вооружение. [1]

По западной классификации ракета получила обозначение SS-N-27 Sizzler (от «sizzle» — шипящий звук, издаваемый кипящим на сковороде маслом). В России и за рубежом (по данным различных СМИ, справочников серии Jane’s и пр.) она обозначалась как Klub, «Бирюза» (Biryuza) и «Альфа» (Alpha или Alfa).

· Назначение

Ракетная система «Клаб» предназначена для поражения надводных кораблей и подводных лодок противника всех типов при ведении боевых действий в условиях сильного радиоэлектронного и огневого противодействия.

· Состав

В состав ракетной системы «Клаб» (Club) входят ракетные (ударные ракетные) комплексы (РК) «Клаб-Н» (Club-N) и «Клаб-С» (Club-S), которые устанавливаются на надводных кораблях и подводных лодках соответственно в качестве ударного ракетного оружия.

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

Ракеты системы во многом унифицированы между собой, но, в зависимости от предназначения и базирования, имеют различные наименования и некоторые отличия. [1]

· Ракета 3М-14

Двухступенчатая крылатая ракета для поражения наземных (береговых) целей подводного (ЗМ-14Э) и надводного (3М-14ТЭ) базирования по внешнему виду, аэродинамической схеме, габаритным характеристикам и двигательной установке аналогична ПКР ЗМ-54Э1 и имеет сходство со стратегической КР ракетного комплекса РК-55 «Гранат» (дальность стрельбы до 3000 км). Отличается фугасной (вместо проникающей) БЧ, подрыв которой осуществляется в воздухе для нанесения максимального ущерба объекту и активной радиолокационной головкой самонаведения АРГС-14Э (ОАО «Радар ММС», г. Санкт-Петербург) с высокоэффективной системой наведения ракеты на цель на конечном участке траектории полета. По этим показателям она превосходит зарубежные аналоги, в т. ч. и американскую Tomahawk, которой может быть поставлена помеха в системе спутниковой навигации GPS. При стартовой массе 2000 кг (БЧ 450 кг) и скорости полета до 240 м/с способна поражать цели на дальности до 300 км. Впервые была показана в февраля 2004 г. на 3 Международной выставке сухопутных и военно-морских вооружений «Defexpo India» (г. Дели). При ее разработке в качестве прототипа была использована стратегическая крылатая ракета «Гранат» (код НАТО SS-N-21 Sampson), предназначенной для вооружения атомных подводных лодок проекта 971, 945, 671РТМ, 667АТ и др. [1]

Таблица 1.1. Основные характеристики 3М54Э1/ТЭ1

· Пусковые установки

Пусковыми установками для пуска всех ракет комплекса «Клаб-С» (Club-S), независимо от их назначения, являются штатные 533-мм торпедные аппараты, а ракет комплекса «Клаб-Н» (Club-N) — подпалубные унифицированные ВПУ (типа ЗС-14Э, «КБСМ», г. Санкт Петербург) модульной конструкции с индивидуальными бронированными крышками (для строящихся) или наклонные ПУ (для модернизируемых) кораблей.

· Модификации

«Клаб-М» (Club-M) наземного базирования разработан для поражения наземных и надводных целей ракетами 3М-14КЭ и 3М-54КЭ соответственно. Дальность обнаружения целей: в пассивном (активном) режиме достигает 450 (250) км, дальность стрельбы по надводным (наземным) целям — до 220 (275) км.

В состав комплекса входят самоходные пусковые установки, ракеты для поражения наземных и надводных целей, машины связи и управления, машины технического обслуживания, комплекс наземного технологического оборудования, учебно-тренировочные средства. Пусковая установка была представлена на МАКС-2007. [1]

2. РАСЧЕТ ДИНАМИКИ ВЫХОДА РАЕТЫ ИЗ ШАХТНОЙ ПУСКОВОЙ УСТАНОВКИ (АНАЛИТИЧЕСКИЙ РАСЧЕТ)

2.1 Исходные данные

Для расчета были приняты следующие исходные данные:

Стартовый вес — 1225 кг;

Длина со стартовым ускорителем — 6,4 м;

Диаметр корпуса — 0,53 м;

Тяга разгонного двигателя — 13 800Н;

Давление в камере сгорания — 7 Па;

Диаметр критического сечения сопла — 0,039 м;

Диаметр выходного сечения сопла — 0,15 м;

Давление атмосферное — 99 600Па.

Так как данные о топливе ракеты 3М-14 отсутствуют, то было выбрано топливо на основе перхлората аммония (68%), полиуретана (17%) и алюминия (15%).

Характеристики топлива:

Удельный импульс — 2460Н;

Газовая постоянная — 290 Дж/(кг*К);

Показатель адиабаты — 1,16;

Температура горения — 3300К. [2]

2.2 Постановка задачи

На основе схемы, изображенной на рис. 2.1, разработать и определить схему расчета параметров выхода ракеты из ШПУ. Для подтверждения качества математической модели и отработки расчетного алгоритма определить один из геометрических параметров ШПУ, а именно ширину отвода газа из «подракетного» пространства и сделать выводы о реальности данного параметра.

Рисунок 2.1. Схема ШПУ с ракетой

2.3 Проектирование ШПУ

Для проектирования шахтной пусковой установки (ШПУ), требуется определить допустимые проектные параметры, которые лимитируются максимальными перегрузками на приборы системы управления и прочностью твердотопливного заряда ракеты с РДТТ.

Примем максимальную перегрузку с равной 15 [2], тогда

.

Таким образом, с учетом этих параметров, определим максимально допустимое значение давление в «подракетном» пространстве :

— Площадь критического сечения сопла

— Массовый расход РДТТ (примем на старте)

кг/с;

где.

— Площадь донного среза

.

— Скорость истечения газа из двигателя

;

где — давление в «подракетном» пространстве.

— Общая тяга ракеты с учетом давления в «подракетном» пространстве

.

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

;

отсюда, допустимое давление .

Составленная система уравнений и последующие уравнения решены в математическом пакете MathCAD.

Теперь, рассмотрим случай, когда объем «подракетного» пространства замкнут, т. е. выход газа из объема отсутствует (рис. 2.2).

Рисунок 2.2. «Замкнутая» схема

— Начальный объем «подракетного» пространства

.

— Теплоемкость газа примем

;

— Температура газа в «подракетном» пространстве

где =+[2].

— Давление будем определять из условия

[4];

где, т.к. система замкнутая — примем, V — скорость выхода ракеты из ШПУ.

— Из уравнения суммы всех сил действующих на ракету, определим V

;

— Соответственно, высоту вертикального выхода ракеты из ШПУ будем определять из уравнения

;

— Составим систему уравнений

Примем .

Решая данную систему уравнений, построим график изменения давления от времени (рис. 2.3) и определим максимальное значение давления.

Рисунок 2.3. Диаграмма

Остальные графики изменения параметров выхода ракеты от времени выглядят следующим образом

Рисунок 2.4. Диаграмма

Рисунок 2.5. Диаграмма

Рисунок 2.6. Диаграмма

Рисунок 2.7. Диаграмма

Рисунок 2.8. Диаграмма

Рисунок 2.9. Диаграмма

Из графика видно, что при t=0.0135секунд будет достигнуто максимальное давление в ШПУ .

Из расчета было определено, что, поэтому необходимо снизить давление за счет установки дополнительного выхода газа из «подракетного» пространства (рис. 2.10).

Рисунок 2.10. Схема ШПУ с отводом газа

Проведем расчет при t=0.0135с, добавив в уравнение

массовый расход отвода газа из «подракетного» пространства. Получим

.

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

Определив допустимый массовый расход отвода, спроектируем оптимальное сечение отвода.

Расчет будем проводить также, при t=0.0135с.

Исходные данные для расчета, полученные при нахождении :

все остальные параметры указаны в общих исходных данных.

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

— Определим среднюю плотность газа из уравнения Менделеева-Клайперона

.

— Динамический коэффициент вязкости из формулы Милликена

[5];

— Кинематический коэффициент вязкости

.

— Площадь выходного сечения отвода газа

.

— Для определения числа Re, найдем осредненную по всему объему скорость газа

.

— Число Рейнольдса Re

.

— Из закона сопротивления Блазиуса, найдем параметр трения

.

— Формула Пуайзеля для труб с некруглым поперечным сечением

[6]

где — коэффициент сопротивления, отнесенный к гидравлическому диаметру; U — смоченный периметр.

Решая систему уравнений (1) — (5), получим, что d₂= 0.586м, т. е. зазор равен 8 мм.

Выводы:

Необходимо отметить, что расчет был проведен в программе MathCAD в интегральной форме с помощью блока Given-Find. Выявлено, что расчет данной модели достаточно неустойчив, т.к. число сложных нелинейных уравнений в системе стремится к числу, при котором MathCAD может давать сбои (5−7 уравнений). Поэтому, необходим подбор параметров расчетной схемы, таких как точность (TOL) и, в особенности, начальные приближения искомых величин. При расчете удалось добиться небольшой погрешности из-за вычислительной схемы Д_max=0.05% (рис. 2.11).

Для сравнения был проведен расчет в дифференциальной форме используя метод Рунге-Кутта, который, как и предполагалось, показал еще большую неустойчивость расчета и, как результат, большую погрешность (Д_max=23%) (рис. 2.12). Так как в системе решаются сложные дифференциальные уравнения с быстро изменяющимися и осцилирующими функциями, то здесь проявилась сложность с подбором расчетного шага. Поэтому, данная расчетная схема в экстремумах функции дала максимальные погрешности.

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

Рисунок 2.11. Диаграмма погрешности вычислений в интегральной модели

Рисунок 2.12. Диаграмма погрешности вычислений в дифференциальной модели

Примечание: Погрешность вычислений Д, определялась из уравнения давлений

— для интегральной схемы

;

— для дифференциальной схемы

.

Из рис. 2.12 видно, что при 0,02 секунды погрешность максимальна (23%), что говорит о том, что параметр p (t) в данной точке не отражает действительность (рис. 2.13), а значит, говорит о не пригодности дифференциальной схемы в данном случае.

Рисунок 2.13. Дианрамма, рассчитанная в дифференциальной схеме

2.4 Расчет параметров выхода ракеты при заданной схеме

Теперь, определим параметры выхода ракеты для схемы с отводом (рис. 2.10) с учетом вычисленных требований по давлению в «подракетном» пространстве.

Составим систему уравнений движения ракеты в ШПУ с учетом отвода газа.

;

На основе полученных данных при решении системы уравнений, построим графики зависимости параметров p, P, V, T₂, U_ист., от времени. Численные результаты приведены в таблице 2.1.

Рисунок 2.14. График p (t)

Рисунок 2.15. График V (t)

Рисунок 2.16. График P (t)

Рисунок 2.17. График T₂(t)

Рисунок 2.18. График (t)

Рисунок 2.19. График H (t)

Рисунок 2.20. График U_ист(t)

Рисунок 2.21. График n_x(t)

Из рис. 2.14 — рис. 2.21 видно, что процессы, протекающие в подракетном пространстве после достижения максимума давления, меняются менее значительно, чем в случае с замкнутым пространством.

Таблица 2.1. Сводная таблица параметров старта

Выводы

Составлена приблизительная схема расчета динамики выхода ракеты из ШПУ с элементами проектирования пусковой установки. В результате расчета определены параметры газа, находящегося в «подракетном» пространстве и параметры движения ракеты от начала движения до ее полного выхода из ШПУ. Кроме того, определены некоторые геометрические характеристики пусковой установки.

Анализ полученных результатов и заключение

В первом приближении расчетное время выхода ракеты из ШПУ составило примерно 0.23 секунды. Для более точного решения требуется уточнить характеристики стартового РДТТ, включить в расчет действие аэродинамических сил при выходе ракеты из шахты. Но следует отметить, что при более тщательном рассмотрении вопроса увеличится число уравнений в системе и, следовательно, сложность вычислений.

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

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

3. МЕТОДИКА ЧИСЛЕННОГО ЭКСПЕРИМЕНТА

3.1 Постановка задачи

Определить и уточнить методику численного эксперимента для решения поставленных задач, в том числе: динамика выхода ракеты из ШПУ и определения АДХ ракеты.

Для реализации численного эксперимента — подобрать оптимальный программный прикладной пакет.

3.2 Общие сведения о проведении численного эксперимента газовой динамики

Первым этапом проведения эксперимента является возникновение или постановка физической задачи. Правильная постановка задачи в значительной мере предопределяет успех всей проводимой работы. Моделирование опыта можно рассматривать как второй этап. Это может быть либо численное моделирование, либо моделирование физических процессов. Каждый из этих методов имеет свои особенности. Численное моделирование дает количественное выражение закономерностей, присущих математической модели. Напротив, с помощью моделирования физических процессов наблюдается аналог действительности или даже сама действительность. Таким образом, полезность расчета ограниченна обоснованностью математической модели.

Численное моделирование задач газовой динамики осуществляется по следующему алгоритму:

1. Построение геометрии модели (тело + расчетная область).

Построение геометрии осуществляется в специализированных пакетах таких, как Solid Works, KOMPAS, Unigraphics, ProEngineer и др. Построенная объемная модель, специального формата, импортируется в сеточный генератор.

2. Построение блочной сетки.

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

3. Подготовка модели (задание свойств среды, начальных условий, граничных условий, параметров решателя).

4. Запуск решателя.

Здесь необходимо проверить «сходимость» решения.

5. Просмотр результатов в модуле Post.

Численное решение, производимое в решателе, представлено на рис. 3.1.

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

На практике дискретизация производных по времени осуществляется почти исключительно с использованием разностных методов. При дискретизации пространственных производных используется, как правило, метод конечных разностей, конечных элементов, конечных объемов [8, 9, 10, 11].

Рисунок 3.1. Структура численного решения

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

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

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

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

Если учесть, что для большинства задач газовой динамики определяющие уравнения являются нелинейными, процесс построения численного решения обычно ведется посредством итераций (используются методы Ньютона, многосеточные методы, метод сопряженных градиентов). Иначе говоря, решение для каждого искомого переменного в каждой узловой точке последовательно корректируется посредством обращения к дискретизованным уравнениям [9, 14, 8, 10, 15, 13].

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

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

· Основные уравнения газовой динамики

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

В наиболее общем случае для задачи газовой динамики требуется решить систему из четырех независимых уравнений, которая носит название системы уравнений Навье-Стокса:

1. Уравнение неразрывности (сохранения массы)

.

2. Уравнение количества движения (сохранения импульса)

где

— тензор напряжений, записываемый в виде

;

— дельта-функция Кронекера

.

3. Уравнение энергии (сохранения энергии)

где

.

4. Уравнение состояния

Для записи соотношений — использованы следующие обозначения: — давление; - плотность; - скорость; - температура; - время; - полная энтальпия; - статическая энтальпия; _ источниковый член для импульса; - источниковый член для энергии; - коэффициент динамической вязкости; - коэффициент теплопроводности; - оператор Гамильтона (набла); - обозначает векторную величину.

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

При рассмотрении некоторых основных дифференциальных уравнений гидродинамики -, можно сделать вывод, что основные переменные подчиняются обобщенному закону сохранения [4, 8]. Если обозначить зависимую переменную, то обобщенное дифференциальное уравнение можно записать в следующем виде:

где

— коэффициент диффузии;

— источниковый член.

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

Анализируя обобщенное дифференциальное уравнение сохранения и саму систему Навье-Стокса, записанную для наиболее общего случая трехмерного нестационарного движения вязкой жидкости, можно видеть, что среди данных выражений присутствуют дифференциальные уравнения в частных производных как первого, так и второго порядка. Дополнительный важный аспект — наличие нелинейной зависимости членов уравнений от переменных.

При историческом развитии динамики жидкости в рассмотрение был введен ряд классов течений, описываемых значительно более простыми системами, чем указанная выше -. Эти различные классы возникают при пренебрежении или ограничении некоторых свойств течений. Для течений, представляющих практический интерес, соответствующая классификация приведена в таблице 2.1. Классификация проведена по двум параметрам — вязкости и плотности. Несжимаемые течения, как правило, ассоциируются с течениями, скорость которых мала по сравнению со скоростью звука (). Наоборот, для сжимаемых течений (, либо разница температур в потоке велика) требуется рассмотреть полное уравнение неразрывности и учитывать полное уравнение энергии.

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

Таблице 2.1. Классификация течений

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

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

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

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

3.3 Компьютерные пакеты для численного решения задач газовой динамики

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

Естественным шагом в эволюции численного моделирования динамики жидкости и газа стало создание расчетных пакетов (CFD-пакетов или комплексов), ориентированных на широкую аудиторию пользователей — научных работников, студентов, инженеров и т. д. В таком виде математический аппарат, заключенный в численные методы, стал действительно универсальным, а с учетом стремительного развития вычислительной техники и мощным средством в проведении численных расчетов по газовой динамике. Кроме того, при использовании CFD-пакетов становится необязательным обладать глубокими знаниями по численным методам, способам дискретизации и т. п.

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

1. «Тяжелые» — комплексы высокого класса, подходящие как для научных, так и инженерных расчетов, способные решать самые сложные задачи с учетом большого количества эффектов и использованием широкого набора математических подходов, в том числе специфических. К классу «тяжелых» относятся лидеры среди коммерческих CFD-пакетов — ANSYS CFX (ANSYS, Inc.), Star-CD (CD-adapco), FLUENT (ANSYS, Inc. совместно с Fluent, Inc.). Все они содержат большое число моделей турбулентности, способны решать задачи различной сложности с учетом горения, химических реакций, многофазных потоков, поддерживают различные типы сеток и т. д.

2. Среднего класса. Предназначены, главным образом, для расчетов инженерного уровня сложности. Набор используемых моделей также может быть достаточно широким. К этому разряду можно отнести COSMOSFloWorks (Solid Works Co.), STAR-CCM+ (CD-adapco), ANSYS FLOTRAN (ANSYS, Inc.).

3. «Легкие» — CFD-комплексы, использующие алгоритмы невысокой точности (используются, например, в качестве учебно-методических), либо имеющие узкую направленность расчета (специально созданные под определенную проблематику).

Следует отметить тот факт, что подавляющее большинство CFD-кодов, реализованных в программах, основано на использовании МКО в различных вариациях.

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

3.4 Общие сведения о методе контрольного объема

Выбранный для примера программный комплекс Ansys CFX v.10 основан на конечно-объемном методе (МКО) решения уравнений гидродинамики. Понимание метода позволит безошибочно провести численный эксперимент, несмотря на то, что в CFX МКО дополнен и утонен рядом других уравнений. Поэтому, можно рассмотреть классический метод контрольного объема, описанный Патанкаром.

В МКО расчетную область разбивают на N-е число непересекающихся контрольных объемов (рисунок 3.2) таким образом, что каждая узловая точка содержится в одном контрольном объеме. Дифференциальное уравнение интегрируют по каждому контрольному объему. Для вычисления интегралов используют кусочно-непрерывные функции, которые описывают изменение зависимой переменной (например, одной из составляющих скорости) между сеточными узлами. В результате находят дискретный аналог дифференциального уравнения. Дискретный аналог представляет собой алгебраическое уравнение,

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

1. Выбор начального приближения или оценка значений Ф во всех узловых точках.

2. Расчет предварительных значений коэффициентов в дискретном аналоге на основе начального профиля Ф.

3. Решение линейной системы алгебраических уравнений, дающее новые значения Ф.

4. Возврат ко второму этапу и повторение процесса до тех пор, пока дальнейшие приближения (итерации) перестанут давать сколько-нибудь существенные изменения в значениях Ф.

Рис. 3.2. Контрольный объём (заштрихованная область), для двумерного случая и обозначения для алгоритма в декартовой сетке

Одним из важных свойств МКО является то, что в нем заложено точное интегральное сохранение таких величин, как масса, количество движения и энергия на любой группе контрольных объемов, следовательно, и на всей расчетной области. Это свойство проявляется при любом числе узловых точек. Таким образом, даже решение на грубой сетке удовлетворяет точным интегральным балансам.

Ansys CFX позволяет проводить расчеты на смешанных сетках, состоящих из различных типов элементов: тетраэдров, призм, клиновидных элементов и гексаэдров (рис. 3.3).

Рисунок 3.3. Ячейки а) — гекса, б) — тетра, в) — призма, г) — пирамида

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

Рисунок 3.4. Схема алгоритма

Выводы

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

Для проведения численного эксперимента был выбран наиболее оптимальный программный прикладной пакет — Ansys CFX, основанный на методе контрольного объема.

4. РАСЧЕТ ДИНАМИКИ ВЫХОДА РАЕТЫ ИЗ ШАХТНОЙ ПУСКОВОЙ УСТАНОВКИ (ЧИСЛЕННЫЙ ЭКСПЕРИМЕНТ)

4.1 Постановка задачи

На основе исходных данных и расчетах, проведенных в разделе 2, разработать стратегию численного эксперимента газодинамической задачи выхода ракеты из шахты. На основе проведенной работы сделать выводы.

4.2 Стратегия постановки численного эксперимета

Газодинамическая задача является временной, поэтому решение ее будет проводиться в зависимости от времени. Для этого необходимо построить несколько моделей, в которых будет меняться объем «подракетного» пространства, в зависимости от движения ракеты. При достижении давления, которое обеспечивает движение ракеты на величину H (величина Н выбирается из условия точности расчета), решается новая модель, начальные условия в которой соответствуют параметрам в последний момент предыдущего расчета. При таком подходе определяются только параметры газа, а движение ракеты определяется из уравнения

и .

· Построение твердотельной модели

Как уже было сказано, модель можно построить в любом графическом редакторе, в том числе KOMPAS, Solid Works, Unigraphics, ProEngineer. В данном случае выбран пакет Solid Works, в котором была создана модель и расчетная область. Так как модель симметрична, то можно ограничиться построением сектора (рис. 4.1).

Конечным результатом построения должны стать сохраненные файлы формата Parasolid (текстовые).

Рисунок 4.1. Схемы твердотельных моделей (нижняя часть)

· Построение сетки

Одним из самых сложных этапов моделирования, является построение сетки расчетной области. Качество сетки влияет на сходимость решаемой задачи и самое главное, на правдивость полученных результатов. Необходимо помнить, что сетка должна иметь больших различий по соотношениям объемов ячеек. При уменьшении объема ячеек — увеличивается точность, но при этом и увеличивается время расчета. Но, зная, что наиболее важной характеристикой точности является величина ячеек на границе с телом, то объем ячеек сетки будем уменьшать от краев расчетной области к границам рассматриваемого тела.

Сетка для схемы изображенной на рис. 4.1 создавалась в специализированном сеточном генераторе ICEM CFD компании ANSYS. Для быстроты построения сетки были использованы тетраэдрические ячейки (рис. 4.2). Но, несмотря на быстроту построения этот метод обладает плохой сходимостью при расчете и снижает скорость самого расчета.

Рисунок 4.2. Тетраэдрическая сетка

· Подготовка решателя расчетного комплекса и проведение расчета

Для экспортированной в Ansys CFX сеточной модели были определены граничные условия. Определение границ изображено на рис. 4.3.

Рисунок 4.3. Определение границ расчетной области

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

1. Модель газа: модель газа, приближенного к параметрам реального воздуха при температуре 25 градусов — Air at 25C.

2. Модель турбулентности: SST с решением полного уравнения энергии.

3. Модель стенки: адиабатическая стенка (Adiabatic Wall) с учетом поверхностных напряжений трения (без проскальзывания) — No Slip.

4. Схема решения уравнений: Временной.

5. Параметры сходимости решения: максимальное число итераций — 30.

6. Временной шаг — 0,005 секунд

7. Граничные условия задавались следующим образом:

8. Граница входа Inlet: статическая температура (Static Temperature) — 2500 К, массовый расход через вход — 5,36 кг/с

9. Свободный выход Open: статическое давление на выходе (Static Pressure for Entrainment) — 99 600 Па, статическая температура (Static Temperature) — 298 К.

10. Граница симметрии — Symmetry.

4.3 Анализ полученных результатов

В результате расчета были получены значения параметров газа в «подракетном» пространстве. Картины параметров газа для t=0.02 c представлены на рис. 4.4. и рис. 4.5.

Рисунок 4.4. Распределение давлений в «подракетном» пространстве

Рисунок 4.5. Линии тока газа в «подракетном» пространстве

Надо сказать, что результаты, полученные при проведении численного эксперимента близки по значениям результатам аналитического расчета на первых секундах (таблица 4.1). Далее, с увеличением времени, увеличивается разность значений параметров газа, полученных при различных расчетах.

Таблица 4.1. Результаты расчетов

Такие неудовлетворительные результаты могли возникнуть из-за ряда причин, в том числе:

1. использовалась модель воздуха;

2. скорость истечения (340м/с) отличается от расчетной (2500м/с);

3. количество уточнений (30) на каждом промежутке времени недостаточно для используемой сеточной модели.

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

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

Выводы

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

5. ОПРЕДЕЛЕНИЕ АЭРОДИНАМИЧЕСКИХ НАГРУЗОК

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

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

Таблица 5.1. Характер и значения скорости ветра

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

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

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

5.1 Постановка задачи

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

В связи с тем, что данная задача мало освещена в литературе и зачастую не поддается аналитическому решению, провести расчет с помощью современного расчетного средства, применяемого для решения задач газовой динамики.

Провести расчет для ракеты выходящей из ШПУ, используя параметры выхода, описанные в первом разделе и ветровой нагрузки.

По результатам расчета, требуется построить эпюры моментов и перерезывающих сил.

5.2 Подготовка и проведение численного эксперимента

· Построение твердотельной модели

Расчетная модель была построена в пакете Solid Works. В данном случае, для исключения влияния стенок на результат, необходимо учитывать площадь заполнения телом расчетной области, требуется, чтоб это параметр был примерно равен 1/8 площади расчетной области.

· Построение блочной сетки

Блочная сетка, для схемы изображенной на рис. 5.1, создавалась в специализированном сеточном генераторе ICEM CFD компании ANSYS.

Рисунок 5.1. Схема моделирования

В результате построения сетки с применением блока «O-grid», т. е. грани ячеек внутреннего блока расположены перпендикулярно к касательным генерируемого тела, была создана расчетная модель (рис. 5.2)

Рисунок 5.2. Сгенерированная блочная сетка

Необходимо отметить, что при разных углах атаки приходится изменять количество ячеек на поверхности тела, так как, к примеру, при изменении угла атаки с 15 на 90 градусов под особое внимание будет попадать сетка на цилиндрической части, а не на носовой. То есть при отработке модели требуется неоднократное корректирование сгенерированной сетки для конкретного расчетного случая.

· Подготовка решателя расчетного комплекса и проведение расчета

Для примера примем, что на ракету в перпендикулярном направлении дует «свежий» ветер, эта характеристика соответствует 10 м/с. Скорость и высоту выхода ракеты возьмем из предыдущего расчета. Для данной схемы (рис. 5.1) проведем расчеты в П.П.П. Ansys CFX.

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

;

.

Для предварительной оценки проведем шесть расчетов с параметрами, указанными в таблице 5.2. Параметры расчетной среды приведены в таблице 5.3.

Таблица5.2. Расчетные параметры

Таблица 5.3. Параметры расчетной среды

Данный расчет проведен с использованием SST-модели турбулентности, так как данная модель более точна и надежна для широкого класса потоков (т.е., потоков подверженных градиентам давления, обтекание профилей), чем стандартная k-w модель турбулентности. [16]

Напомним, что жидкость, находящаяся в расчетной области, при M << 1 будет считаться несжимаемой.

· Получение результатов

В результате расчета были получены теневые картины обтекания ракеты. Приведем такие картины для расчетного примера № 6 (рис. 5.3, рис. 5.4).