Расчет нестационарной гидродинамической нагрузки на крышку пускового устройства при подводном выбросе ла избыточным давлением газа Текст научной статьи по специальности «Физика»

УДК 532.5.031;629.7

Расчет нестационарной гидродинамической нагрузки на крышку пускового устройства при подводном выбросе ЛА избыточным давлением газа

© А.В. Плюснин1,2

1МГТУ им. Н.Э. Баумана, Москва, 105005, Россия 2 ОАО «ВПК «НПО машиностроения», Московская область, г. Реутов, 143966, Россия

Предложен приближенный способ расчета давления и силы на крышку пускового устройства (ПУ) при подводном выбросе ЛА, основанный на методе граничных элементов (МГЭ).

Ключевые слова: подводный выброс ЛА, нестационарные гидродинамические силы, метод граничных элементов.

Введение. Одним из способов подводного выброса ЛА является использование работы расширения газов, находящихся в ПУ под давлением, превышающим наружное гидростатическое давление [1-2]. При выходе кормы ЛА из ПУ газы оттуда вырываются в окружающую среду и формируют в области между кормой ЛА и срезом ПУ газовый пузырь (ГП), сначала расширяющийся, затем схлопываю-щийся, причем избыточное давление газов передается окружающей воде. Этот кратковременный процесс сопровождается интенсивной знакопеременной нагрузкой на близлежащие конструкции. В данной работе рассматривается силовое воздействие на крышку ПУ, находящуюся в открытом положении (рис. 1).

Рис. 1. Воздействие избыточного давления на крышку ПУ

Предлагаемый подход основан на МГЭ [3-7] и является более экономным, чем использование методов конечных разностей или конечных элементов. В то же время данный подход позволяет провести существенно более полный анализ задачи, чем оценки, опирающиеся на фундаментальное решение уравнения Лапласа в пространстве (задача Рэлея о расширении сферической полости в жидкости) и теорию присоединенных масс [8].

Формулировка приближенного подхода к решению задачи на основе МГЭ. Характер временных процессов в рассматриваемой физической задаче не требует учета сжимаемости воды, а также вязких и вихревых эффектов. Поэтому движение жидкости будем описывать на основе математической модели потенциальных течений однородной несжимаемой жидкости. Это дает следующие соотношения в области, занятой жидкостью [8]:

div V = 0,

V = grad ф,

АФ = 0, (1)

дф V2 р рн

— + — + — + gz = ¿-н-,

д 1 2 Рж Рж

где V = V, V - вектор скорости частиц жидкости; ф - потенциал ско-

, д2 д2 д2 рости; А = —2 н--2 ^--2 - оператор Лапласа в декартовых коор-

д х д у д г

динатах х, у, г; рж - плотность жидкости; g - ускорение силы тяжести; р - абсолютное давление в жидкости; рн = ратм +рж gH -гидростатическое давление на уровне верхнего среза ПУ; ратм - атмосферное давление; Н - глубина погружения ПУ; 1 - время.

В формуле (1) перейдем от абсолютного давления к избыточному давлению 5р = р - рн + рж gz, что дает

5р = -Р,

V2 ^

чд1+2 У

(2)

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

сти, при этом вдали от ПУ ф^ 0, V ^ 0, Ър ^ 0. С другой стороны, считаем, что поверхность конструкции, в которой находится ПУ, достаточно велика. Заменим эту поверхность неограниченным плоским экраном, на поверхности которого должно выполняться условие непротекания:

о п

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

Чтобы избежать в последующем интегрирования по поверхности, уходящей в бесконечность, применим зеркальное отображение области, заполненной жидкостью, относительно плоскости г = 0. Для этого отобразим также корпус ЛА, присоединенный к его корме ГП и саму открытую крышку ПУ (рис. 2).

♦ г

Рис. 2. Область решения задачи

Движение ЛА при выходе из ПУ считается заданным, т. е. путь Ь , скорость V и ускорение А известны как функции времени I. Естественно предполагать, что присутствие крышки ПУ слабо влияет на характер расширения и схлопывания ГП. Это предположение, возможно, будет нарушено при очень больших давлениях газа в момент

выхода ЛА из ПУ, однако на практике стараются таких ситуаций не допускать.

В рамках сделанного предположения будем считать, что задача динамики ГП, являющаяся намного более трудной, каким-либо образом решена, например, используя соответствующие экспериментальные данные или простые численно-аналитические оценки. Итак, считаем, что известен закон движения границ ГП во времени и давление среды внутри ГП. Располагая этими данными, в области жидкости между крышкой ПУ, поверхностью ЛА и ГП (с учетом зеркального отображения) требуется решить для ряда моментов времени уравнение Лапласа относительно потенциала скорости ф с соответствующими граничными условиями.

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

д 1

поверхности. Поскольку для решения уравнения Лапласа известны представления в виде граничных интегральных уравнений [9], то подходящим методом для решения задачи является МГЭ [2-4]. Част-

дф

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

д 1

функция. В данной работе для нахождения дф ставится отдельная

д 1

краевая задача. Дискретизация этой задачи по методу МГЭ приводит к СЛАУ с той же матрицей (и другой правой частью).

Вывод граничного интегрального уравнения (ГИУ). Охватим

поверхности ЛА ЕЛА, ГП ЕГП , крышки Е^ и их отражения Е ЛА, ЕГП и Е*- воображаемой сферой Е^ радиуса Я, с центром в начале координат (рис. 2). На крышке ПУ зафиксируем регулярную точку поверхности М (хм, ум, гм) и окружим ее воображаемым шаром Ое

достаточно малого радиуса в. Пусть £Е - часть поверхности шара, находящаяся в жидкости и Е^ = Е|с^р)\ Ов. Введем в рассмотрение также подвижную точку Р (х, у, 2) и обозначим расстояние между

точками М и Р через м =^(х- Хм)2 +(у - Ум)2 +(- 2м)2 .

В области ,в, ограниченной поверхностями ЕЛА = ЕЛА и Ела,

ЕГп =ЕГППиЕГП, Екр = ЕкриЕЕе , Ефункции ф, ^ф и —

д 1 гмр

являются гармоническими. Применяя формулу Гаусса - Остроградского, выведем вторую формулу Грина для этих функций [9]:

0= | (фа

Дш ;е

1 1 л^

---аф

ГМР ГМР у

а □=

Д« ;е

Ф§гаё-

' МР у

-

1

V ГМР

§гаё Ф

а □=

1 0Ф

ГМР 0 п

(р )-Ф(Р )±

( л V

V ГМР у

(4)

где Ф = ф, Оф; 2^ ;е=2 ла и 2 гп и 2^ у 2<? у 2еу 2д

1

На поверхности 2к при функции Ф и - имеют по-

рядок О

С 1 л

V

интегралы

ОФ 0 а - и —

0 п 0 п

С 1 л

V гмр у

имеют порядок О

С 1 л

'МР

С 1 л

Дш

. Тогда

V1« у

Г 0ф(р)аъ и Г Ф(Р)-°- аъ стремятся к

г-~ 0 п 2Д 0 п V гмр у

2д Гмр 0 п

С 1 л

нулю в пределе ^ ш .

Н 2 110 На поверхности 2е имеем -= — и —

Гмр в 0 п

= —12 . Следовательно, потенциал простого слоя |

V гмр у

= §гаё-

1 МР

ГМР ГМР

1 0Ф

2 гмр 0 п

( р ) аъ

стремится к нулю при в^ 0, а для потенциала двойного слоя, учитывая регулярность точки М, имеем

гф(р )а

) v ;0п

С 1 л

с18 =

V гмр у

1Ф(М )£С л *+ДФ(Р )-ф(М )]0

0 п V гмр у 2

С 1 л

с18 =

V гмр у

= -2 лФ(М ) + 0 (1).

Переходя в формуле (4) к пределам и в ^ 0, получим

основную интегральную формулу Грина:

I

ш

ц

2 лФ(М)+ j

ула и угп и укр

1 ЭФ

д

'MP

(P )-Ф(р)— д n W V ' д n

1

V rMP

dS = 0. (5)

В этой формуле оба интеграла J

1 дФ

у( +) rMP дn

^кр

(P) dS

и

г 3 Г 1 ^

I Ф(Р)— - dS мажорируются абсолютно сходящимися не-

у(+) 3 п V ГМР )

^кр

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

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

з Г

"37

Отсюда имеем

д n ч

J

у rMP д П

кр

= 0.

(6)

(P) dS = 0.

Вычисление потенциалов простого слоя (ППС) и двойного слоя (ПДС) в осесимметричном случае. Поверхности ЛА и ГП в каждый момент времени t считаются известными. Будем полагать их формы осесимметричными, что в основном и бывает на практике. В цилиндрических координатах r ,9 , z ( х = r cos 9 , y = r sin 9 ) эти поверхности удобно описывать параметрически как функции длины

дуги ^ : r = Гла (St),

z = z

ла

(s, t)

и r = r

гп

(s, t), z = Zrn (s, t). Тогда

дФ,

входящие в правую часть ГИУ (5) Ф( Р) и —— (Р) также могут быть

представлены как функции только от ^ и I. В этом случае интегрирование ППС и ПДС по поверхностям ЛА и ГП можно представить как интегрирование по дуге меридионального сечения границы и по углу 9 [2-4]. Интегралы по 9 сводятся к эллиптическим интегралам, для которых известны простые аппроксимации [10].

Таким образом, процедура вычислений будет следующей. Меридиональные сечения границ ЛА и ГП разбиваются на криволинейные , , 3Ф/ ч

ГЭ. Функции Ф(Р) и —— (Р) аппроксимируются по значениям в

центрах элементов. Интегрирование по дуге ядер ППС и ПДС проводится с помощью, например, квадратурных формул Гаусса, исходя из следующих соотношений, учитывающих зеркальное отражение относительно плоскости 2 = 0 (суммирование по индексу ] означает суммирование по всем криволинейным ГЭ, составляющим границу ЛА и ГП в области 2 > 0 ):

-1- Г

9 тт Л

1 0Ф

2 — 2 гп и 2 ла ГМР 0п т ;+1

(Р) ёЪ =

( з+У2)I КПспС(гм , гм, г г и) ^

2 —

Г ф-0.

0 п

2гп и 2ла

с 1 л

с18 =

V гмр у

+1

ЕФ(] +1/2 ) I 4дсм)гм , 2М, г 2 (<*) г'(°) *'(<*))^

где

КППС (гм , 2м,г,2) = КППС (гм , 2м,г,2) + КППС (гм , 2м,г, 2), К (сим)( г ггг 'г ') =

= КПДС (гм , 2М,Г,г',2') + КПДС (гМ , 2М , г,~2,—,2' , к ППС (гм , 2м, г, 2 ) = 2 га ( гм , 2м, г, 2 )к(к 2), КПДС (гм , 2м,Г,г',2') =

= {в (гм , 2м , г , 2, г ', )(к 2 )-2 '[к(к2 )-Е(2 )]}А ( , 2м , Г, 2),

1

а (гм , 2м, г, 2 ):

—

В (

гм , 2м, г, ^ г , 2

к2 =-

у1(гм + г)2 + (м - 2)2

' л 2 г [ 2 ' ( Гм - г )-г ' ( 2м - 2 )]

)=

(гм - г)2 + (м - 2)2

4 гмг

(гм + г) +(м -

2

С.+1/2 " V 2

к( )= I

0

V 2

С . + а.+1

d ф

- к2 мп2 ф

Е(к2)= | А/1-к281п2фdф .

Вычисление потенциалов двойного слоя (ПДС) в пространственном случае. В ГИУ (5) остается вычислить только интегралы

3

вида IФ(Р^

Г 1 ^

-'кр

V гм.р )

dS. Учтем наличие в задаче еще одной плос-

кости симметрии: у = 0. Выполним аппроксимацию той части поверхности крышки, которая находится в октанте х>0, у<0, г>0, набором плоских треугольных элементов и зеркально продолжим это разбиение относительно плоскостей симметрии г = 0 и у = 0 . Искомые функции Ф( Р) аппроксимируем значениями в центрах ГЭ. Таким образом, задача свелась к вычислению ПДС на плоских треугольных элементах.

Пусть точки , ^2, ^з задают вершины треугольного ГЭ и его ориентацию: из конца вектора нормали п (положительная сторона ГЭ) последовательный обход вершин Ау, А2, А3 должен быть виден

А

3 п

против часовой стрелки. Интеграл I —

* Ям

а1а2 а3

Г 1 ^

dS равен телесно-

V ГМР )

му углу, под которым из точки М видна положительная сторона ГЭ [9], т. е. знак этого телесного угла противоположен знаку скалярного произведения МА1 • п. По формулам сферической геометрии [11] получим

3п

I -

* Г11

а1а2 а3

Г 1 ^

V ГМР )

dS = а(тт-р1-р2-Рз),

где

а = sign

ign (МА1 • п);

МА

а =

2

аi = arccos (а. • ак );

* 1

Sin

0,5 (а. + ак-а; )

sin аг

Рг = 2arcsin(Ъ. Ък)

и индексы г, . и к образуют циклическую перестановку из индексов 1, 2 и 3.

Составление СЛАУ, порядок вычислений. Искомыми величинами ГИУ (5) относительно Ф = ф являются значения потенциала

ф(М) в центрах плоских треугольных элементов, покрывающих часть поверхности крышки в октанте х>0, у<0, г>0. Коэффициенты при этих значениях образуют матрицу СЛАУ. Вектор правых ча-

3

стей СЛАУ вычисляется в результате интегрирования ф и- вдоль

3п

границ меридионального сечения поверхностей ЛА и ГП как объяснено выше. Для численного решения СЛАУ использовался стандартный прямой метод [12].

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

3

Основное значение имеет решение ГИУ (5) относительно Ф =-.

3 ^

Искомыми величинами являются значения -^ф- (М) в центрах треугольных ГУ. Из ГИУ (5) следует, что матрица СЛАУ получается той же, что и в случае Ф = ф . Вектор правых частей второй СЛАУ вычис-

Г3ф>

вдоль границ

3ф 3 Г3ф

ляется в результате интегрирования - и — -

31 3 п V 31

меридионального сечения поверхностей ЛА и ГП.

тг 0ф

После нахождения значений v и- вдоль поверхности крышки

0г

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

^ (г) = - | Ъp (Р, г) п сК.

-кр

Распределение абсолютного давления на крышке определяется формулой

Р = ЪР+Рн-Рж ^

Рассмотрим для примера процесс, протекающий по схеме (см. рис. 1). Зададим параметры движения ла (Ь (г) V (г), А( г)) и избыточным

давлением бргп (г) в ГП согласно рис. 3. Разбиение поверхности крышки на треугольные элементы показано на рис. 4, где изображены центры ГЭ в проекциях на координатные плоскости Х2 и ху .

Рис. 3. Параметры ЛА при выходе из ПУ и давление в ГП

Процесс расширения и разрушения ГП описать математически довольно сложно, поскольку имеют место эффекты взаимопроник-

Г% 4

1

г £ --* №№ (*№ № Чх» *«*

>1 Л Г* л Е* 5>

4 < • е- о + о ■ » -1 " -а 4

'> <5 * ^ <> * ■ г О' + & 3 *

Л9-*-* Ьш Ф ■ ФФ > * % о4 г!

с"5»,» -о $ ф <

о * 4Щ < * * О > « > «

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

Рис. 4. Разбиение поверхности крышки на треугольные ГЭ

В процессе со сферической (цилиндрической) симметрией влияние жидкости падает обратно четвертой (второй) степени расстояния,

тогда как влияние •ддф падает лишь обратно пропорционально рас-

д ^

стоянию (корню квадратному от расстояния). Следовательно, в рассматриваемой задаче главное значение имеют размеры ГП и величи-

д

на- на его границе (и близлежащей поверхности корпуса ЛА). Вед ^

дф ттг

личина - на поверхности 111 находится непосредственно из соот-

д ^

ношения (2). Остальные граничные условия на ГП легко получить (в рассматриваемом примере), пользуясь решением задачи Рэлея [9]. На боковой поверхности ЛА, считая ее цилиндрической, имеют место условия (3), (6), потенциал ф можно взять из решения задачи обтекания тела, производную -дф- вычислить дифференцированием (учи-

д ^

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

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

Результат расчета силы (по нормали к базовой плоскости крышки) показан на рис. 5.

Рис. 5. Расчетное значение силы на крышку ПУ

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

ЛИТЕРАТУРА

[1] Апальков Ю.В., Мант Д.И., Мант С.Д. Отечественные баллистические ракеты морского базирования и их носители. Санкт-Петербург, Галея Принт, 2006, 216 с.

[2] Дегтярь В.Г., Пегов В.И. Гидродинамика подводного старта ракет. Москва, Машиностроение, 2009, 448 с.

[3] Власова Е.А., Зарубин В.С., Кувыркин Г.Н. Приближенные методы математической физики. Москва, Изд-во МГТУ им. Н.Э. Баумана, 2001, 700 с.

[4] Banerjee P.K., Butterfield R. Boundary Element Methods in Engineering Science. London, McGraw-Hill Book Company (UK) Limited, 1981, 494 p.

[5] Brebbia C.A., Telles J.C.F., Wrobel L.C. Boundary Element Techniques. Berlin, Springer-Verlag, 1984, 523 p.

[6] Boundary Integral Equation Method: Computational Applications in Applied Mechanics. Ed. by Cruse T.A., Rizzo F.J. New York, ASME, 1975, 390 p.

[7] Chuang Y.K., Ehrich O. On the Integral Technique for Spherical Growth Problems. Int. J. Heat and Mass Transfer, 1974, vol. 17, рp. 945-953.

[8] Лойцянский Л.Г. Механика жидкости и газа. Москва, Наука, 1978, 736 с.

[9] Соболев С.Л. Уравнения математической физики. Москва, Наука, 1966, 444 с.

[10] Abramowicz M., Stegun I.A. Handbook of mathematical function. New York, Dover, 1965, 832 p.

[11] Бронштейн И.Н, Семендяев К.А. Справочник по математике для инженеров и учащихся втузов. Москва, ГИТТЛ, 1957, 608 с.

[12] Demmel J.W. Applied Numerical Linear Algebra. Philadelphia, SIAM, 1997, 430 p.

Статья поступила в редакцию 27.06.2013

Ссылку на эту статью просим оформлять следующим образом:

Плюснин А.В. Расчет нестационарной гидродинамической нагрузки на крышку пускового устройства при подводном выбросе ЛА избыточным давлением газа. Инженерный журнал: наука и инновации, 2013, вып. 7. URL: http://engj ournal.ru/catalog/mathmodel/aero/847.html

Плюснин Андрей Владимирович окончил Днепропетровский государственный университет по специальности «Гидроаэродинамика» и аспирантуру кафедры «Волновая и газовая динамика» МГУ им. М.В. Ломоносова. Канд. физ.-мат. наук, доцент кафедры «Вычислительная математика и математическая физика» МГТУ им. Н.Э. Баумана, заместитель начальника отела ОАО «ВПК «НПО машиностроения». Автор ряда работ по нестационарным задачам гидрогазодинамики и теории упругости. е-mail: andrey.plusnin@gmail.com