Переходные характеристики массивного тонкого крыла в ограниченном потоке Текст научной статьи по специальности «Механика и машиностроение»

Вычислительные технологии

Том 16, № 2, 2011

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

Е. П. Лукащик Кубанский государственный университет, Краснодар, Россия e-mail: lep_9091@mail.ru

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

Ключевые слова: переходная характеристика, задача Вагнера, удлинение крыла, экран, дискретные вихри, дозвуковой поток.

Введение

Определение нестационарных гидродинамических сил на движущихся телах необходимо при изучении динамики и устойчивости тел и решении задач оптимизации и управления движением. Движение крыльев или крыльевых систем вблизи границы произвольной формы представляет собой сложный неустановившийся процесс. Течение вокруг крыла в общем случае имеет трехмерную структуру. В математической физике хорошо известен подход к исследованию общих нестационарных процессов, основанный на использовании так называемых переходных характеристик. Переходная характеристика динамической системы определяет зависимость от времени параметров системы при переходе от одного стационарного состояния в другое. Такой переход можно моделировать, задавая системе возмущение входного сигнала в виде "ступеньки" или единичной функции Хевисайда. В теории систем автоматического регулирования [1] обычно с помощью преобразования Лапласа по времени строится передаточная функция Ш(в), по которой затем получают переходную характеристику с помощью обратного преобразования Лапласа

h (t) = L

-i

W (s)

(t).

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

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

1. Движение массивного тонкого крыла конечного размаха в несжимаемой жидкости вблизи экрана

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

где Н (¿) — функция Хевиеайда, уо — амплитуда кинематических возмущений, — скорость набегающего потока;

— отсутствие перепада давления в вихревом следе за крылом

со следующими краевыми условиями: — непротекание крыла

(1)

и на боковых кромках крыла;

— непротекание твердого экрана

ду

0, г = -к, Ух, Уу,

дг

и нулевыми начальными условиями

<Р (х, У, г, 0) = (х, у, г, 0) = 0.

В тексте используются обозначения: Бр : {|х| < а, |у| < Ь} — проекция прямоугольного крыла на плоскость ХОУ, : {|х| > а, |у| < Ь} — проекция вихревого следа па плоскость ХОУ, а — полухорда, Ь — полуразмах крыла.

Краевое условие (1) включает параметры движения крыла как твердого тела: / (¿) — смещение центра масс, а (¿) — угол поворота крыла относительно оси, перпендикулярной скорости потока и проходящей через центр масс в горизонтальной плоскости, в (¿) _ угол крена.

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

= Р +

= МУ~ МУ° -

= мх- Мх0 - М°х (2)

с нулевыми начальными условиями

/(0 ) = |(0) = а(0) = ^(0)=/!(0) = |(0) = 0.

Здесь М* — масса крыла; Зу и ,1Х — моменты инерции относительно осей ОУ и ОХ; Р — главный вектор гидродинамических сил реакции; Му — момент тангажа, Мх — момент крена; Я0 — главный вектор внешних сил динамичеекого возмущения, Му0, Мх0 — главные моменты внешних сил относительно осей ОУ и ОХ, М°, МХ — внешние моментные пары. Для этих величин справедливы следующие выражения:

Р = 11 (Р- - Р+) Яо = \\ Яо (^ пП)

МУ = С (Р- - Р+) , Мх = V (Р- - Р+) ,

Му0 = !I СЯо (^ пП) Мх0 = Ц 'ПЯо (^ пП) (3)

где я0 (х, у) — поверхностное напряжение внешних динамических сил.

Для построения базового уравнения динамики массивного крыла над экраном воспользуемся методом гидродинамических особенностей [5], согласно которому воздействие на возмущенный поток жидкости тонкого крыла можно заменить воздействием системы п-образных вихрей, непрерывно распределенных по поверхности крыла и моделирующих в каждой ее точке необходимый скачок давления. Такие вихри называют присоединенными.

Элементарный п-вихрь состоит из двух бесконечно длинных вихревых нитей единичной интенсивности, но разного знака, проходящих параллельно оси ОХ на бесконечно малом удалении друг от друга и соединяющихся между собой бесконечно малым отрезком вихревой нити также единичной интенсивности [6].

Потенциал элементарного п-вихря с вершиной в точке (£,п,я) определяется как интеграл потенциалов диполей, расположенных вдоль бесконечной нити

те

4п { у/(х - Л)2 + (у - п)2 + (г - я)2 (1+ (х-0

4тг (г-я)2 + (у~л)2 \ ^ _ ^ + {у _ ^ + (г _

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

РпН (х - - п,г) = Рп (x,í,У, п, г, 0) - Рп п, г, -Щ.

Предположим, что нестационарность движения возникла с момента £ = 0, Интенсивность присоединенных п-вихрей 7 (£,п,£)> моделирующих крыло, в силу нестационарного характера обтекания крыла меняется со временем. Изменение интенсивности по схеме Бирнбаума—Кюсснера [4] сопровождается возникновением свободных вихрей, вместе с частицами жидкости перемещающихся вниз по потоку и образующих вихревой след за крылом.

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

Р(х,У,г,£) = // (1(^,П,£)РпЬ(х - С,У - П,г)-

г

/д7

V, т)ржН(х -£-Ь + т,у-г1, г)<1т )(1з. (4)

о

Потенциал в форме (4) удовлетворяет уравнению Лапласа и условиям в точках следа и экрана.

Величина скачка давления на поверхности крыла определяется интенсивностью присоединенных вихрей

7 (£,п,£)

1 др др уте д£ дх

V- ~Р+

рУте '

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

Р (£) = рУте [[ 7 (е,п,£) ¿8,

My (t) = pvх п, t)

Mx (t) = pv0

nj (C,n,t) ds.

(5)

Уравнение для вычисления интенсивности присоединенных вихрей получим из условия непротекания на крыле (1), используя представления для потенциала возмущенного течения вблизи экрана в виде (4) и гидродинамических характеристик в виде (5) без учета влияния внешних сил:

Y (C,V,t) R (x - C,V - n,h) -

dj

(С,П,Т) R (x - С - t + T,y - n,h) +

dr\ds = v0H

где R (x,y,h) — функция влияния вихря для плоской твердой границы

R(x,y,h) = {x,y,z)

(6)

dz

z=о

Интегральное уравнение (6) для интенсивности присоединенных вихрей 7 (С,п,£) является гиперсингулярным по переменной вдоль размаха. Значение интеграла с такой особенностью понимают в обобщенном смысле Адамара—Манглера [6]:

J (y - п) £^0 \J (У - П)

a У a y+£

7 (v) , 27 (у)

( \2 1 <г

(У - п) £

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

Понижения особенности подынтегральной функции в уравнении (6) можно добиться, применяя интегрирование по частям вдоль размаха крыла и учитывая, что 7 (х, у, £) на боковых кромках в любой момент времени принимает нулевое значение, В результате получим второй вариант уравнения для определения интенсивности присоединенных вихрей 7 (х, у, £)

dj

т-(£,ri,t)k(x-£,y-r],h) -

d2j д'цдг

, п,т) k (x - С - vo (t - т) ,y - n,h) +

dj

1 yy

м* j„ 2 JT

dr | dr? = —VQH ^ 00

(x,y) e Sp, (7)

где

k (x, y, h)

1(1/ | y^TF

4n | y \ x

y

4h2 + y2

t

t

ху (8к2 + х2 + у2)

у/ж2 + у2 + Ш2 [4/12 (4/г2 + ж2 + у2) + х2у2}

Учитывая связь вихревой интенсивности с перепадом гидродинамического дни. ю-

д7 (х, у, г)

пия, определяем, что искомая функция д(х,у,1) = -^—'— конечна на выходной

ду

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

Интегральное уравнение (7) примем в качестве базового уравнения динамики массивного тонкого крыла вблизи экрана.

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

/ (г)

а (г) =

ру0

М*

(г - т) ^ / / 7 (С,п,т) ^

рУо

Jy

рУо

М*

рУо

{г-т)(1т / / - 1

а - т) ¿т Л г,Л

а - т) &т

дг)

в (г)

РУо

Jx

(г - т) dW / VI (С,пп,т) ^р

рУо

(г - т) dт п

>.д1 (С,п,т)

дц <18р-

(8)

Для прямоугольного крыла верпы следующие выражения (к0 — малая толщина,

ро

т 7 4а3Ь , 4аЬ3 М = р0П0 ■ Аао, ■1у = р0П0-^-, .1х = р0П0

3

Введем безразмерные величины

_7 _х_уЛЬ-аЬ

7=—, х = —, у = т, А = -, ¿ =—, га

Уо а Ь а Уо

и запишем уравнение (7) в безразмерном виде

РоЬр ра

3

I Ь - Уо 1г= ь0 = — а Уо

д7 д'ц

(С,П,г) к (х - С,\ (у - п) ,Ь) -

<927

д'цдт

(С,п,т) к (х - С - (г - т) (у - п) }к) +

dт\ dСdn = -у0н(г).

(9)

г

г

г

г

В уравнении (9) черточки над безразмерными величинами опущены, В дальнейшем

Б Г -1 < х < 1,

Бр Л -1 <у< 1.

Коэффициенты подобия на основе безразмерного решения уравнения (9) определяются по следующим формулам:

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

р—ь Зр зр

— коэффициент момента тангажа

коэффициент момента крена

Мх 1 г г 1 г г 2ду&г1,г)

^ = в УУ77 —^—'

где Б — площадь крыла,

— относительные перемещения центра масс крыла

г г

1а=И {<~Т)ЛТ II ~<«'Т) ^1^11

0 Зр о Зр

— угол тангажа

г г

± I (, - г) Лг И {7 К, „, г) А, I (( - т) Лг II

о Зр о

— угол крена

г г

нг)=^1{г-т)ат11гп*>>т)**р1{г~тит IIг]2д1%п,т)^-(1о)

о Зр о Зр

В результате дискретизации нестационарного процесса во времени нестационарную задачу (9) сведем к ряду рекуррентно связанных стационарных задач [8]:

9 (С,'П,гт)

Аг / 3уп

К (х, у, г/, к, А, г, г) — -—г] ( 1 + Зх£ +

4т V 2

dСdn = -уоН (г)-

К (х, у, г/, к, А, г, /с) — -—г] М + + ^ '

// ]>> (С,п,г*)

к=1

(и)

где

R (x,£,У,п,Ь,Л,r, к) = = к (х - £ - (гг - гк), Л (у - п), Ь) - к (х - £ - (гг - 4-1), Л (у - п), ь).

Для каждого расчетного момента времени 1Г двумерное интегральное уравнение (11) решается численно. Для построения численной схемы используется классический метод дискретных вихрей [5]. Поверхность крыла представляется в виде прямоугольных элементов, для этого размах крыла разбиваем равномерно на 2М частей, а хорду на 2М частей. При численном решении интегральных уравнений с ядром Коши для класса решений, неограниченных на одном конце интервала и конечных на другом, эффективной является следующая дискретная схема расположения присоединенных вихрей (£3, п^) и контрольных точек (хг,у^), в которых выполняется условие непротекания:

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

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

По изложенной схеме рассчитаны переходные характеристики массивного крыла конечного удлинения в ограниченном потоке. Для значительных величин относительной массы крыла удлинения Л = 4 на больших расстояниях от границы и при Ь = 0.2 расчеты хорошо согласуются с результатами работы [9]. Учет массы крыла приводит

а б

Рис. 1. Влияние высоты полета массивного крыла (а) и массы крыла вблизи границы (б) на подъемную силу

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

2. Нестационарное обтекание массивного крыла бесконечного удлинения ограниченным потоком сжимаемой жидкости

Ь

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

(1 М2) 9^ 9^ - 1 9^ 2М 9^ 12

^ ' дх2 дг2 с2 <9£2 с дхдV

V»

где с — скорость звука, М =--число Маха набегающего потока,

с

Краевые условия на профиле, твердом экране, вихревом следе и начальные условия идентичны условиям задачи для крыла конечного размаха в ограниченном потоке несжимаемой жидкости, рассмотренной в разделе 1,

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

(й2 ~ ' ^ (й2 ~ у-В безразмерном виде (за единицу длины выбрана полухорда профиля а, за единицу скорости — V») уравнение (12) примет вид

+ 0 = + (13)

Условие непротекания в точках профиля следующее:

др ч df ¿а г, -

ог ¿ъ ¿ъ

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

[р- (x, ъ) - р+ (x, ъ)]

7 ОМ) =-•

pv»

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

(£, п)

определяется выражением

Ро (х,£,г,п,ъ) =

^т.

1 ( д7 (С,т) - Л)у(г - т) -М (х - С - г + т) - лУ

= / о arctg-*-^-

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

Рн (х, С, z, г) = (x, С, z, 0, г) - (x, С, z,-2h, г).

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

+1 г

/ Л 1 Г „ Г&У&т) -1 о

arctg ■

(г - т) - м2 [(х - С - г + ту + z2

г2 + (х-0 (х-с-г + т)

— arctg ■

^ + 2К)^(г - т)2 - м2 [(х - С - г + т)2 + ^ + 2К)2 {г + 2 к)2 + (х-£)(х-£-Ь + т)

dт,

подстановка которого в условие непротекания профиля приводит к интегро-дифференциальному уравнению

+1 г г +1

^ /<£/ = + ( ¿т /7(£,г)(1 + (14)

-1 о

о -1

где

к (х, г) = к1 (х, г) - к2 (х, г) , к1 (х, г) =

к2 (х, г) =

ф2-м2 (х-гу х (х - г)

(х (х - г) + 4Н2) (г2 - м2 [(х - г)2 + 4Н2] - 4М2н2) - 8Н2 (г2 - м2 [(х - г)2 + 4Н2])

[(х (х - г) + 4Н2)2 + 4Н2 (г2 - м2 [(х - г)2 + 4Н2])] ф2 - м2 [(х - г)2 + 4Н2]

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

Су (*) = /7 (е, г) сту у) = - / ет (е, *) с

-1

г 1

-1

г 1

/(^ = ¿1 Ц-г) 17 I ЫМ^т.

о -1 о -1

Для численного решения уравнения (5) проведем дискретизацию нестационарного

г

1

1

Рис. 2. Зависимость коэффициента подъемной силы (а) и коэффициента момента тангажа (б) от безразмерного времени; Н = 0.2, М = 0.5

изменение кинематических параметров, на равные интервалы АЪ. Непрерывное изменение интенсивности 7 (х,Ъ) заменим ступенчатым. В силу такой замены образование свободных вихрей будет происходить также в дискретные моменты времени. Производ-д^у (х, Ъ)

ную-г—^— аппроксимируем конечной разностью

дъ

<97 (х, и) _ 7 (х, и) - 7 (х, ¿¿-1) дЬ А1 '

Учитывая нулевое значение 7 (х,Ъ) при Ъ = 0, такими преобразованиями интегро-дифференциальное уравнение (14) сведем к системе рекуррентных интегральных уравнений

+1

/ 7 (С,и) \к (х -i.tr- гг-1) - —Аг ■ (1 + Зж01 =2тгу0Ни)+ II т

+1

-1

-1 г—1

1=1

УЧ(£,*г) к(х-£,и-и)-к(х-£,и-11_ 1) + -А1-(1 + ЗхО (15)

т

Численное решение системы (15) выполнялось на основе метода дискретных вихрей [5]. Выбор узлов £3 и точек коллокации хг определялся поведением искомой функции х, Ъ) на кромках профиля. Учитывая связь интенсивности присоединенных вихрей с перепадом давления, имеем, что искомая функция 7(х,Ъ) является неограниченной на передней кромке и конечной на выходной кромке крыла. Для данного класса решений эффективной будет следующая дискретная схема:

На рис. 2 показано влияние массы крыла на поведение во времени коэффициента подъемной силы Су и коэффициента момента тангажа Сту {К{= 0.5, Ь = 0.2).

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

Заключение

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

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

[1] Попов Е.П. Теория линейных систем автоматического регулирования и управления. М.: Наука, 1989. 304 с.

[2] Ворович И.И., Бавешко В.А., Пряхина О.Д. Динамика массивных тел и резонансные явления в деформируемых средах. М.: Научный мир, 1999. 247 с.

[3] Бавешко В.А., Глушков Е.В., Зинченко Ж.Ф. Динамика неоднородных линейно-упругих сред. М.: Наука, 1989. 345 с.

[4] Некрасов А.И. Теория крыла в нестационарном потоке. М.: Изд-во АН СССР, 1947. 258 с.

[5] Белоцерковский С.М., Скрипач Б.К., Табачников В.Г. Крыло в нестационарном потоке газа. М.: Наука, 1971. 768 с.

[6] Эшли X., Лэндал М. Аэродинамика крыльев и корпусов летательных аппаратов. М.: Машиностроение, 1969. 318 с.

[7] Гахов Ф.Д. Краевые задачи. М.: Наука, 1972. 640 с.

[8] Лукащик Е.П. Апериодические движения крыльев над плоской твердой границей // Самолетостроение. Техника воздушного флота. Харьков: Выща школа, 1988. С. 14-20.

[9] Гур-Мильнер С.И. Численный метод определения аэродинамических производных движущегося вблизи земли крыла произвольной формы в плане // Труды ЛКИ. 1977. Вып. 115. С. 41-46.

[10] Белоцерковский С.М., Кочетков Ю.А., Томшин В.К. Динамика движения тела при учете нестационарности обтекания //Механика твердого тела. 1974. № 3. С. 35-43.

[11] Ефре мов И.И., Липован П.С., Лукащик Е.П. Переходная характеристика подъемной силы при апериодическом обтекании тонкого профиля ограниченным сжимаемым потоком // Экологический вестник научных центров ЧЭС. 2008. № 4. С. 51-57.

Поступила в редакцию 25 мая 2010 г., с доработки — Ц октября 2010 г.