CC BY
106
Челябинский физико-математический журнал. 2019. Т. 4, вып. 1. С. 65-75.
УДК 621.455; 629.76.085.5 Б01: 10.24411/2500-0101-2019-14106
ЧИСЛЕННОЕ ИССЛЕДОВАНИЕ ГИДРОДИНАМИКИ ПОДВОДНОГО СТАРТА РАКЕТ
В. И. Пегов
Южно-Уральский научный центр УрО РАН, Миасс, Россия
Государственный ракетный центр им. академика В. П. Макеева, Миасс, Россия
Предложены математическая модель и метод расчёта гидродинамики старта ракет из затопленной шахты подводной лодки. С помощью конформного преобразования неосесимметричная область течения в шахте переводится в осесимметричную, и трёхмерная краевая задача сводится к ряду осесимметричных задач, решаемых методом Галёркина в сочетании с методом конечных элементов. Нестационарное взаимодействие заполняющей шахту жидкости и погружённой в неё ракеты выражается через зависимости присоединённой массы от геометрических параметров ракеты и шахты: ширины кольцевого зазора, формы ракеты, положения границ раздела сред «газ-жидкость».
Предложенный алгоритм расчёта присоединённой массы реализован в виде программы на ЭВМ. По результатам верификации численного моделирования с данными аналитических решений подтверждена достаточная достоверность разработанной математической модели гидродинамики ракеты в ограниченном объёме жидкости.
Ключевые слова: ракета, старт, гидродинамика, присоединённая масса, конформное преобразование, метод конечных элементов.
Старт ракет с подводных лодок является одним из наиболее сложных видов старта, так как в реальных условиях, создаваемых ходом подводной лодки, волнением моря и высоким давлением окружающей среды, ракеты испытывают повышенные нагрузки, значительно превосходящие нагрузки на воздушной траектории. Имеющиеся экспериментальные и расчётные данные, а также опыт эксплуатации спроектированных в АО «ГРЦ Макеева» ракет морского базирования позволяют полагать, что наиболее безопасным и эффективным является способ старта, осуществляемый запуском маршевого жидкостного двигателя ракеты в затопленной шахте, который не требует применения других корабельных стартовых устройств, обеспечивает наименьшие стартовые осевые и боковые нагрузки на ракету и сохраняет возможность управления движением ракеты в воде. В литературе данные по гидродинамическим силам, действующим на ракету при выходе её из затопленной шахты, представлены весьма ограниченно, и изучение их с помощью современных численных методов остаётся актуальным.
Форму ракеты обычно представляют в виде удлинённого тела вращения, состоящего из кругового цилиндра со сферическим притуплением носка и плоским кормовым срезом, шахту — в виде круглой трубы с глухим днищем. Движение ракеты происходит вдоль прямой, соединяющей оси ракеты и шахты, со скоростью и0. При решении задачи будем считать жидкость идеальной и несжимаемой, а вызванное течение жидкости в шахте — потенциальным.
Очевидно, что гидродинамическая сила, действующая на ракету, удовлетворяет соотношению
d
Р = -~т (т*щ), (1)
dt
где t — время, т* — присоединённая масса ракеты, определяемая из выражения для кинетической энергии жидкости Т:
2
Т = . (2)
В отличие от движения тела в безграничной жидкости при движении в ограниченном пространстве величина производной ^т1 не равна нулю из-за влияния стенок. Учитывая это, из соотношения (1) получаем выражение для гидродинамической силы Р
dm* duo
Р = - т*1л. (3)
Величина присоединённой массы т* является функцией расстояния Ь между осями ракеты и шахты т* = f (Ь), которая меняется с течением времени. Очевидно, что
dm* dm* dЬ
dt dЬ dt dt
Из соотношения (3) с учётом последних равенств получаем
2 dm* duo
F = - m*lf - (4)
Таким образом, в соответствии с (4) для определения гидродинамической силы F необходимо найти величину присоединённой массы и закон её изменения m* = f (b) в зависимости от расстояния между осями ракеты и шахты. Введём декартову х,у, z и цилиндрическую т,в, z (x = r cos 9, y = r sin 9) системы координат, ось z совпадает с осью симметрии ракеты, начало координат помещается на днище шахты, ось x направим вдоль прямой, соединяющей оси ракеты и шахты (рис. 1). Форму ракеты зададим зависимостью Ri(z) текущего радиуса образующей ракеты от координаты z, через Яц обозначим радиус цилиндра ракеты, радиус шахты — через R2, через Si, S2 — поверхности ракеты и шахты, S3, S4 — поверхности торца ракеты и днища шахты. Положительными будем считать внешние к поверхностям нормали.
Очевидно, что потенциал течения ^(r,z,9) удовлетворяет уравнению Лапласа в цилиндрической системе координат
1 д Íтдф\ + 1 dV + = о (5)
т дт \ дт J т2 дв2 dz2
и следующим граничным условиям:
д^ = — u0 cos в cos y на S1} (6)
дп
где y — угол между касательной к образующей ракеты и её осью,
д^ = 0 на S2, S3, S4. дп
Выражение для кинетической энергии жидкости (2) запишем в виде
S
М2 + í д^2 + (2
дх J \ду) \дz
dT =2 ОД, (7)
Рис. 1. Физическая область течения г* = х + 1у
где — интеграл Дирихле, т — занимаемый жидкостью объём.
Из соотношений (2) и (7) следует выражение для присоединённой массы
т* = -2Б(ф) = —2 [[ <^-¿8.
и0 и0 Л дп
я
(8)
Учитывая свойство инвариантности интеграла Дирихле при конформных преобразованиях, целесообразно производить вычисления присоединённой массы т* с помощью конформного преобразования в области более простого вида, например, в осесимметричной области, когда будут совпадать оси ракеты и шахты. В плоскости течения введём комплексную переменную г* = х + %у.
С помощью дробно-линейного преобразования [1]
ш = —а2-
г* — а1 'г* — а2
(9)
неконцентрические окружности поверхностей $1 и $2 перейдут в концентрические в плоскости ш окружности поверхностей $1 и $2 (рис. 2).
В отображении (9) а1 и а2 — координаты точек, симметричных одновременно относительно окружностей $1 и Б2:
а1
(а — /Ж—Г)
Я1
а2
(л + /Ж—Г) Я1, а
. (ГО)
При этом отрезок положительной оси х между окружностями $1 и $2 переходит в отрезок положительной оси и между окружностями $1 и $2, радиус окружности
Рис. 2. Преобразованная область течения ш = и + ;т
остаётся равным а радиус шахты в преобразованной области течения находится по формуле
а2(а1 + &) (11)
R'
R2
При выводе формулы (11) использовались выражения (10).
Введём в плоскости ш координаты p,fi,z,u = р cos в, v = р sin в. Уравнение Лапласа (5) для преобразованной области течения примет вид
1 5
Р^" +
1 52'
0,
рН2 dp\dpj р2 H¡ дв2 dz2 где коэффициент Ламэ Hi с учётом преобразования (9) равен
й2(й2 - a i)
H
р2 + 2ра2 cos в + a2
(12)
(13)
Из выражения (9) можно установить связь между аргументами 9 и в на поверх-
ностях Si и , где z* = R ieie, ш = R ^,
cos 9
1 + A cos в A + cos в ,
dz *
du
_ VA2—Г
p=R A + cos в'
Граничные условия (6) в плоскости ш запишутся в виде [1]
д' д' дп' дп
dz*
du
P=Ri
^-- (1+ A cos в) '
-u0V A2 - 1—-— cos y на S ,,
(A + cos в)2
^ = 0 на S2, S3, S4. (15)
дп' 2 3 4
Граничное условие (14) можно переписать в виде
= -u0VA2 - 1 cosYf (ß) на S[, (16)
f = dp = öp дн' dp
P=Ri
где
f (ß) 1 + A COS ß
(А + cos в )2'
Для решения краевой задачи (12)—(15) с помощью метода Галёркина [2] запишем
следующее условие, эквивалентное данным соотношениям:
2п
Я / Ят\ Я / 1 Ят\ Я / _ Ят\ 1
¿pdß V dpdz = 0. (17)
z р ^ 0
д f dp\ д /1 ö / 2öp dp + dßVPd^/ + dZ vp 1 dZ
Вариационная задача (17) решается с помощью представления искомого потенциала в виде тригонометрического ряда Фурье
те
р(р, г, в) = -иоЯц ^ Фк(р, г) ^ кв, (18)
k=i
где Фк(р,г) — неизвестные функции ряда Фурье.
Приведённое разложение потенциала в ряд Фурье допустимо вследствие того, что заданная на поверхности ракеты S1 нормальная скорость жидкости представи-ма в виде
f = Y1 fk COs kß3,
k=1
где с учётом формулы (16) коэффициенты ряда Фурье равны
2п
fk = ^ if cos kßdß = (-1)k fcäk-1(1 - ä?)uo cos Y, äi = (19)
n J Ri
0
Ряд Фурье для функции Hf, входящей в соотношение (17), имеет следующий вид:
H = Y + Е hk cos kß, k=1
где коэффициенты разложения равны
2п __k
hfc = 1 У H2 cos kede = ^^ ( -d - ^ (k^l —d2 + l) , k = 0,1, 2,...
° (20) В равенстве (20) введены обозначения
й2 («2 - ai) , 2ра,2 " -, d =
p2 + ä2 p2 + of
Значение к0/2 определяется равенством (20) при к = 0
ко_ 2
(1 - (2)3/2
(1 - а2)2
1+4
«2
Ч)
1 - Г
а
Р< 1. а2
(21)
Выражение (17) с учётом равенств (18)-(20) после интегрирования по частям и введения функционалов
к =
г р
дфкЛ 2 + дфк др ) + 2 V дг
+ к2 Фк р2
р(р(г + Д у Фкр(И (22)
1о
перепишется в виде 8Fk = 0, к = 1, 2, 3,..., т. е. решение трёхмерной краевой задачи сведено к решению ряда двумерных (осесимметричных) задач по минимизации функционалов (22), которую мы осуществим с помощью метода конечных элементов [2; 3]. Для этого область течения в плоскости р,г разобьём на треугольники. Обозначим через Ф^ неизвестные функции, каждая из которых определена на своем г-м треугольном элементе, тогда искомые функции Фк(р,г) можно представить в виде суммы
п
Фк (р,г) = ^ Ф^(р,г), к =1, 2, 3,..., (23)
г=1
где п — общее количество треугольных элементов.
Функции Ф^\р, г) в пределах элементов аппроксимируются через узловые неиз-
вестные {Фк} и функции формы
N
(0
следующим образом [2]:
Ф
(0
N
(0
{Фк} , к = 1, 2, 3,..., г = 1, 2, 3,... ,п.
С учётом выражения (23) функционалы (22) имеют вид Fk = ^ Fk , где
г=1
F,
(0
1 2
А(*)
дФ
(^2 Но( д Ф« \2 к2 (Ф«
др
+ -г
2 \ дг
+
р2
Р<Л(г) + Дк Фкг)р(Ь(г). (24)
Заметим, что ниже в основном используются обозначения и сокращения, принятые в работе [2], причём
В
(0
10
0 ^ 2
в
(0
д^
(0
др дМ
(0
дг
где к0 приведено в равенстве (21);
в
(0
— матрица градиентов функций формы.
В матричном виде соотношение (24) перепишем следующим образом:
F,
(0
1 2
А«
{Фк}
Т
Т
вк1) вк°
в
(0
{Фк} + -2 {Фк} р2
Т
N
(0
Т
N
(0
{Фк}
?(Л(г)+
2
с
2
к
2
к
N
(О
т
{Фк} к = 1, 2, 3,
(25)
Для минимизации функционалов (25) на множестве узловых значений {Фд} необходимо выполнить равенства (п — общее количество элементов в области течения)
дЯд
д {Фк}
д {Фк}
0, к = 1, 2, 3,...
i=1
Последние соотношения представим в матричном виде
ал
£
i=1
К
(0
а {Фк}
где использованы обозначения
т
{Фк } + Р
,(0
0, к = 1, 2, 3,
(26)
К
(0
т
В
В
(0
э<1А(1:) +
т
Л«
Л«
р
Р,
(0
(0
К
(0
Если ввести матрицы [Кк] = ^
i=1
примут вид
[Кк] {Фк} = {Рк}
N
{Р }
т
р
"Е\Р>
i=1
к = 1,2, 3,...
то соотношения (26)
(27)
Таким образом, задача минимизации функционалов (25), а тем самым и краевая задача (12)-(15) окончательно сведена к решению ряда систем линейных уравнений (27) относительно неизвестных узловых значений потенциалов Фд.
В результате с учётом равенств (8), (16), (18), (19) формула для безразмерной присоединённой массы ракеты
т *
приводится к виду
т
рпЯцЬ
т* = -
1 ™ г
X Е Фк —кЯ к=1 {
где
— , Я1
По
Я1 (г)
Я
X Я _ г
ц
Я
ц
ц
Описанный способ вычисления присоединённой массы ракеты был реализован в виде программы для ЭВМ, разработанной А. В. Федотовым. По результатам верификации численного моделирования с аналитическими решениями движения шара в жидкости, заключённой внутри сферической концентрической оболочки [4], и движения эллипсоида вращения в жидкости, заполняющей софокусную эллиптическую полость [5], максимальные отклонения не превысили 3.5 %. Поэтому можно заключить, что численное моделирование достаточно достоверно описывает и течение в ограниченном объёме шахты. Для разбиения преобразованной области течения р, г использовалось до 1800 треугольных элементов, взаимосвязанных в 1001 узловой точке. Результаты расчётов, полученные по разработанной программе, представлены в виде графиков на рис. 3-7.
0,5 Ь
Рис. 3. Зависимость присоединённой массы цилиндра от его смещения в трубе для различных удлинений
Рис. 4. Зависимость присоединённой массы от ширины кольцевого зазора для фиксированных смещений Ь (п = 1.1)
Рис. 5. Зависимость присоединённой массы ракеты от положения свободной нижней границы жидкости
и
Рис. 6. Распределение присоединённой массы в зависимости от вида и положения нижней границы жидкости
Были выполнены расчёты цилиндрического тела с прямыми срезами при поперечном движении его в трубе в зависимости от величины смещения Ь = Ь/ (К2 — Кц) (0 < Ь < 1.0) для различных его удлинений Л (0.5 < Л < 20) П = К2/Кц = 1.1 (рис. 3). Присоединённая масса цилиндра зависит как от удлинения Л, так и от величины смещения Ь. Чем больше удлинение, тем ближе присоединённая масса по величине и характеру изменения к аналитической кривой, полученной для случая цилиндра бесконечного удлинения [1], и наоборот. При приближении цилиндра к стенке трубы (Ь ^ 1.0) присоединённая масса начинает резко расти.
На рис. 4 приведены графики изменения присоединённой массы ракеты в зависимости от безразмерной ширины кольцевого зазора шахты 8 = (К2 — Кц)/Кц. В расчётах было принято: Л = 10, 0.1 < 8 < 1.0, п = 1.1. Присоединённая масса значительно возрастает при уменьшении 8 от 0.3 до 0.1. При заданном значении 8 присоединённая масса тем больше, чем больше величина смещения Ь.
Была проведена оценка влияния перемещающейся вверх по задонному объёму и кольцевому зазору шахты нижней свободной границы жидкости (рис. 5). Из графиков рис. 5 следует, что при перемещении свободной границы в задонном объёме (при 80 < 0) присоединённая масса остаётся практически постоянной; при дальнейшем перемещении свободной границы вверх по кольцевому зазору присоединённая масса ракеты уменьшается приблизительно по линейному закону.
На рис. 6 приведены распределения присоединённой массы по длине ракеты для различных положений свободной границы (80 = 0; 6; 10; 14; 16), а также для случая, когда свободная граница, находящаяся на уровне среза кормы (при 80 = 0), заменяется на твёрдую. Из сравнения распределений присоединённой массы по длине ракеты для двух видов нижней границы жидкости видно, что кривые
т.
15
распределении существенно отличаются друг от друга в кормовой части ракеты
Из графиков рис. 7 следует, что влияние формы носовой части ракеты на присоединённую массу незначительное и тем меньше, чем больше удлинение ракеты.
Результаты выполненных численных исследований позволяют построить зависимости присоединённой массы ракеты от геометрических параметров, которые используются для расчёта действующих на ракету гидродинамических сил при движении ракеты в шахте.
10
7.5
п-/./ -п/Ы WH IfJrl -Ъ-О-сфсра
>.=20
>.=10
0
0.5
Список литературы
Рис. 7. Зависимость присоединённой массы от формы носовой части тела для различных его удлинений
1. Капанкин, Е. Н. Удар кругового цилиндра в жидкости с внешней жёсткой или свободной границей в виде окружности. Учёт влияния границ потока в аэродинамических трубах / Е.Н. Капанкин // Тр. ЦАГИ. - 1969. - Вып. 1154. - 25 с.
2. Сегерлинд, Л. Применение метода конечных элементов / Л. Сегерлинд. — М. : Мир, 1979. — 392 с.
Коннор, Дж. Метод конечных элементов в механике жидкости / Дж. Коннор, К.Бреббиа. — Л. : Судостроение, 1979. — 263 с.
Ламб, Г. Гидродинамика / Г.Ламб. — М.; Л. : ГИТТЛ, 1947. — 928 с. Mikelis, N. E. The added mass of confocal ellipsoid / N. E. Mikelis, A. G. Parkinson, W. G. Price // International Shipbuilding Progress. — 1981. — Vol. 28, no. 318. — P. 40-47.
Поступила в редакцию 19.12.2018 После переработки 16.02.2019
Сведения об авторе
Пегов Валентин Иванович, доктор технических наук, профессор, ведущий научный сотрудник, Южно-Уральский научный центр УрО РАН, Миасс, Россия; главный научный сотрудник, Государственный ракетный центр имени академика В. П. Макеева, Миасс, Россия; e-mail: [email protected].
Chelyabinsk Physical and Mathematical Journal. 2019. Vol. 4, iss. 1. P. 65-75.
DOI: 10.24411/2500-0101-2019-14106
NUMERICAL STUDY OF HYDRODYNAMICS OF UNDERWATER MISSILE LAUNCH
V.I. Pegov
South Ural Scientific Centre of the Ural Branch of the Russian Academy of Sciences, Miass, Russia
Academician V.P. Makeyev State Rocket Centre, Miass, Russia [email protected]
A mathematical model and a calculation method for hydrodynamics of a missile take-off from a submerged silo of a submarine are offered. A non-axisymmetric flow area in the mine is transformed into an axisymmetric one by a conform mapping, and a 3D boundary value problem is reduced to a series of axisymmetric problems that are solved by Galerkin's method with using the method of finite elements. A nonstationary interaction between a liquid filling the silo and a missile loaded in it is expressed through dependencies of apparent mass on the following geometry parameters of the silo and the missile: a width of a circumferential gap, the missile shape, and a location of the gas-liquid interface.
The provided calculation algorithm of the apparent mass is implemented as a computer program. The resulted verification of the numerical simulation with data of analytical solutions confirms sufficient fidelity of the developed mathematical model of the missile hydrodynamics in the restricted liquid volume.
Keywords: missile, take-off, hydrodynamics, apparent mass, conformal transformation, method of finite elements.
References
1. Kapankin E.N. Udar krugovogo tsilindra v zhidkosti s vneshney zhyostkoy ili svobodnoy granitsey v vide okruzhnosti. Uchyot vliyaniya granits potoka v aerodinamicheskikh trubakh [An impact of a round cylinder in a liquid with clear or free external boundary in the form of a circle. Accounting for the influence of flow boundaries in aerodynamic tunnels]. Trudy TsAGI [Proceedings of Central Aerohydrodynamical Institute], 1969, iss. 1154. 25 p. (In Russ.).
2. Segerlind L. Primeneniye metoda konechnykh elementov [Application of the finite elements method]. Moscow, Mir Publ., 1979. 392 p. (In Russ.).
3. KonnorJ.J., BrebbiaC.A. Finite Elements Techniques for Fluid Flow. London, Boston, Newnes — Butterworths, 1976. 309 p.
4. Lamb H. Hydrodynamics. New York, Dover Publ., 1932. 738 p.
5. MikelisN.E., Parkinson A.G., Price W.G. The added mass of confocal ellipsoid.
International Shipbuilding Progress, 1981, vol. 28, no. 318, pp. 40-47.
Accepted article received 19.12.2018 Corrections received 16.02.2019
