__________УЧЕНЫЕ ЗАПИСКИ ЦАГИ
Том XXX 19 99
№1 — 2
УДК 629.735.33.015.3.025.1
О СОКРАЩЕНИИ РАЗМЕРНОСТИ ЗАДАЧИ РАСЧЕТА ПЕРЕХОДНЫХ ФУНКЦИЙ КРЫЛА ПРИ МАЛЫХ ДОЗВУКОВЫХ СКОРОСТЯХ
Н. Н. Глушков
В рамках численного метода дискретных вихрей в линейной постановке рассмотрена задача расчета переходных функций для тонких несущих поверхностей в несжимаемой среде. Показано, что приращения аэродинамических сил вследствие схода вихрей в поток определяются изменением скоростей в окрестности передних кромок. Это обстоятельство позволяет сократить размерность задрчи и свести ее к расчету аэродинамического нагружения крыла в стационарно/I потоке, но с переменной во времени круткой сечений. Приведены примеры расчета переходных функций для крыльев. В двумерном потоке получены точные решения и дано сравнение с теоретическими результатами [1]. При переходе от дискретных распределений вихрей к непрерывным для производных переходных функций получены соответствующие интегральные уравнения.
В нестационарной аэродинамике крыла одна из основных задач — расчет переходных функций, характеризующих процесс установления течения, возникшего при внезапном изменении граничных условий. В течение переходного периода с кромок крыла сходят свободные вихри, которые изменяют его обтекание. Поэтому в каждый расчетный момент времени необходимо определять поле скоростей на поверхности крыла и решать краевую задачу [2]. Покажем, что поле скоростей можно вычислять не на всей поверхности крыла, а только в окрестности передних кромок, что позволяет сократить размерность задачи. Этот результат был получен в работе [3].
1. Дискретные вихревые системы крыла. Пусть в момент времени
7 = 0 тонкое деформированное крыло из состояния покоя внезапно начинает движение с постоянной скоростью У„. С этого момента начинается переходный процесс к стационарному обтеканию. Предположим, что стационарное обтекание крыла моделируется системой дискретных прямоуголь-
ных вихрей, расположенных на базовой плоскости крыла (7 = 0), с равномерным шагом Ах по оси X (рис. 1 ,а). Временной интервал, в течение которого происходит установление течения, разобьем на равномерные отрезки At = Ах/¥Х.
В течение рассматриваемого отрезка времени Д/ на крыле происходит изменение циркуляции присоединенных вихрей, сопровождающееся образованием свободных вихрей, которые сносятся набегающим потоком.
Проведем вертикальную плоскость через контрольные точки в произвольном сечении крыла (например, через сечение А—А на рис. 1 ,а). На
(»)
а)
1 <1 > г
О о
О о о о
о 0 1 я < » X ° 1 0 1
1) Э
о о о
о о о О о
■о о о О О
і)
Рис. 1
временном интервале Л/
(и = 1, 2, ...) в рассматриваемой плоскости (рис. 2,а) будет находиться (п +1) вихревая система, включая систему присоединенных вихрей и п систем свободных вихрей. Для наглядности на рис. 2,а свободные вихри, образовавшиеся в различные промежутки времени, разнесены по высоте, хотя они находятся на одном уровне. Выполнение условий непротекания позволяет определить циркуляцию присоединенных вихрей.
Произведем вычитание вихревых систем, соответствующих двум соседним временным интервалам д/”^ и Д^"-1\ Получим систему диполей, распределенных с кусочно-постоянной плотностью (рис. 2,6). Эта система отвечает за приращение циркуляций присоединенных вихрей на временном интервале А^п\ Плановая проекция дипольных систем крыла показана на рис. 1,6.
2. Развитие дипольных систем в начальный период времени. На
первом временном интервале (п = і) дипольная система крыла находится под воздействием набегающего потока. Условия непротекания в рассматриваемом сечении крыла (индекс сечения опускаем) имеют вид
ДК.(1)+<х,. =0, 1 = 1,2,...,#,
(1)
-АТ}1)
! ^
-дг/-» _ _і Ц З-трС'
с° с° с°
м
В) Система диполей
±Л1}
±4Г/ге"^
$ С0) С0) с°
'УХ® со с°
Рис. 2
а) 1 — присоединенные вихри; 2 — свободные вихри; 3 — контрольные точки; 4 — поле нормальных скоростей от набегающего
потока; 5 — поле нормальных скоростей от вихрей;
б) 1 — присоединенные вихри; 2 — свободные вихри; 3 — кон-
трольные точки
где Д — распределение приращений безразмерных скоростей, вызванных приращениями циркуляций ДгР^; а,- — распределение местных углов
атаки; N — число контрольных точек в сечении.
При и = 2 дипольная система крыла находится под воздействием ранее образовавшегося слоя диполей, смещенного на расстояние Ах = К00 • Дг. Условия непротекания
ДГ/2)+ДР£>=0, / = 1,2,...,#, (2)
с учетом (1) принимают вид
АУ}2) + АК0(1) = О, Д^.(2) - ам =0, I = 2, 3, ..., N. (2а)
Индекс «0» соответствует точке, расположенной перед крылом.
При п = 3 новая дипольная система крыла находится под воздействием двух предыдущих, которые сместились по потоку вместе со своим полем скоростей. Условия непротекания АУ^ + АУ^ + АУ^ = 0 с учетом (2) приобретают вид
Д г/3) + ДР0(2) + Д уМ = 0, ДК.(3) =0, / = 2, 3, ..., N.
Условия непротекания для временного интервала имеют вид
при этом во всех контрольных точках, кроме первой, нормальные скорости равны нулю.
На основании (3) можно утверждать, что, начиная с третьего временного интервала, дипольные системы, образовавшиеся ранее, экранируют крыло, т. е. не пропускают возмущения к его поверхности. Открытыми для возмущений остаются только контрольные точки, прилегающие к передним кромкам. Следовательно, приращения аэродинамических нагрузок в переходный период определяются изменением скоростей вблизи передних кромок.
Распределение скоростей вдоль кромок крыла является функцией одной геометрической переменной, изменяющейся по размаху крыла. Только эта функция входит в граничные условия, и через нее могут быть выражены циркуляции крыла и рассчитаны все аэродинамические характеристики. Приходим к основному результату: размерность задачи расчета переходного процесса может быть сокращена на одну переменную.
3. Эквивалентная форма сечений крыла. Для каждого временного интервала можно построить крыло, которое в стационарном потоке будет иметь такое же распределение циркуляции присоединенных вихрей, как и рассматриваемое крыло в нестационарных условиях, т. е. будет эквивалентно исходному в смысле аэродинамического нагружения.
Возвращаясь к рис. 1, напомним, что стационарное обтекание моделируется системой П-образных вихрей. В то же время, если систему П-образных вихрей сдвинуть вниз по потоку, изменить знак циркуляции и сложить с исходной, получим дипольную систему. Следовательно, скорости в контрольных точках дипольной системы и эквивалентной ей системы П-образных вихрей связаны соотношением
Как видно, местные углы атаки суммируются. Если исходное крыло
в-Гм)«,=*!
или
(4)
где Г0еч — скорость перед кромкой крыла от П-образных вихрей.
І
При и = 1, согласно (1), имеем
У=1
является плоским и начинает движение под углом атаки а, то
(У,- -Уп)^ =-«•/, 7 = 1, 2, Ы, и нормальные скорости в сечении изме-еЧ
няются по линейному закону. Следовательно, форма сечения — парабола. Разность скоростей в окрестности передней кромки ограничена
-а, т. е. течение здесь не имеет особенности.
Таким образом, на первом временном интервале исходное плоское крыло имеет такое же аэродинамическое нагружение, как крыло с параболическими профилями и круткой. Крутка эквивалентного крыла обеспечивает режим плавного (безударного) обтекания передней кромки. При п = 2 с учетом (1), (2) и (2а) получаем
К/2) = АУ}1) + ДГ,(2) = -СЦ - ДГ0(1), ^(2)=ДГ,(1)+ДГ/2)=-а,+<хм, / = 2,3,...,#.
Для эквивалентного крыла с учетом (4) имеем
7=1
поэтому с точностью до константы нормальные скорости соответствуют распределению местных углов атаки. Следовательно, при гг - 2 эквивалентное крыло возвращается к форме исходного крыла и отличается от него только местными углами атаки (круткой) сечений. С течением времени крутка крыла восстанавливается.
Крутка эквивалентного крыла может рассматриваться как искомая функция времени. Через нее определяются циркуляция крыла и аэродинамические коэффициенты.
4. Влияние деформаций крыла. Предыдущий анализ показывает, что деформация сечений крыла участвует в переходном процессе только на первых двух интервалах времени. В дальнейшем решается задача о возмущениях на передних кромках (3) или задача о крутке эквивалентного крыла. Рассмотрим подробнее влияние деформаций.
Решение стационарной задачи можно проводить как по схеме П-образных вихрей (рис. 1 ,а), так и по схеме диполей плюс один П-образный вихрь в конце каждого сечения крыла (рис. 1,в). Последний вихрь содержит всю циркуляцию сечения, поскольку сумма циркуляций дипольных пар равна нулю. Если с помощью крутки устранить циркуляции в сечениях крыла, то на таком крыле в стационарном потоке останутся только диполи.
При п = 1 обтекание крыла также моделируется системой диполей (рис. 1,6), причем граничные условия непротекания в обоих случаях одинаковы. Следовательно, распределения диполей, изображенных на рис. 1,6 и в, должны совпадать. Это возможно только в том случае, если диполи, примыкающие к задней кромке (рис. 1,6), равны нулю. Приращения цирку-
ляций присоединенных вихрей связаны с циркуляцией вихрей в стационарном потоке следующей зависимостью:
АГ/1)=ЕГй.У» *- = 1,2, (5)
М
Здесь Г817 — циркуляции вихрей в стационарном потоке. Свободные вихри имеют циркуляцию обратного знака ДГ^. = -ДГ^, / = 1, 2, ..., N. Равенство нулю циркуляций дипольных пар, прилегающих к задним кромкам
N
крыла, ДГдР = Гз1у = 0 приводит к тому, что при смещении дипольного 7=1
слоя по потоку при п = 2 за крылом не окажется диполей.
На втором временном интервале на крыле под первым слоем диполей возникает второй слой. Поскольку за крылом диполей нет, то для выполнения условий непротекания (2) циркуляции дипольных пар второго и первого слоя должны быть взаимно уничтожены:
ДГ/2)+ДгД=0, / = 1,2, ...,#,
где ДГ^ = 0.
Следовательно, перед крылом не возникнет дополнительных возмущений и переходной процесс заканчивается. Суммарная циркуляция присоединенных вихрей распределена по закону
г/2) = ДГ/!) + ДГ/2) = ДГ/!) - ДгД{, / = 1,2,...,#.
С учетом (5) имеем = Г8и, / = 1, 2, ..., N, причем при п = 2 устанавливается стационарное распределение циркуляции. Свободные вихри, образовавшиеся при п = 1, поглощаются новым свободным вихревым слоем с циркуляцией противоположного знака
4Гй-1+4гй=-[дгй+лгР)] = °. ' = '• 2. N.
где ЛГ"’о“0'
5. Крыло с плоскими сечениями. Задачу расчета переходных функций деформированного крыла (отклонение закрылков, носков, прогибы профилей и т. д.) полезно разделить на две независимые задачи. Одна предполагает расчет переходного процесса для крыла с плоскими сечениями и круткой, обеспечивающей такой же закон распределения циркуляции по размаху, как у исходного крыла. Вторая задача проводится в стационарном потоке для деформированного крыла с той же дополнительной крут-
кой, но обратного знака, когда циркуляция в сечениях равна нулю. Суммирование распределений циркуляций этих крыльев дает решение нестационарной задачи. Уменьшения времени расчета здесь можно добиться за счет уменьшения количества вихрей вдоль хорд по сравнению со случаем деформированных сечений. Кроме того, для плоских сечений крыла решение базовой задачи о возмущениях на передней кромке можно начинать не с третьего временного интервала, а со второго и рассматривать не приращения циркуляций, а текущее состояние дипольной системы крыла. Действительно, пусть а5 — угол атаки рассматриваемого сечения крыла. Тогда, суммируя скорости, представленные в уравнениях (1), (2а), (3), получаем
п /1-1 т п
т=1 т=\к=\ т= 1
/ = 2, 3,..., ІУ, п = 2, 3, ... .
приращение нормальной скорости перед крылом;
полная нормальная скорость перед крылом, учитывающая вклад от набегающего потока;
нормальные скорости от дипольной системы кры-
п
ла с циркуляциями = Елг/т)-
т= 1
у(п)=_у(п-1). у(п)=0} / = 2>3, ...,#, п = 2,3,.... (6)
Для эквивалентного крыла, согласно (4), получим (У\ - У0 = -У^п Вводя обозначение (^о)еч =(*о +а5')еч и учитывая условие непротекания (У[ +а5)еч =0, получаем
(^=^1 п = 2, 3, ... . (6а)
Следовательно, нормальная скорость перед кромкой эквивалентного крыла с запаздыванием на один шаг по времени равна нормальной скорости перед кромкой рассматриваемого крыла.
Введем обозначения:
т
к=1
т=1 п
у(п)
т=1
Тогда вместо (3) имеем
6. Базовая задача переходного процесса. Вводя обозначения номеров сечений и полагая в (6) и (6а) У^п ^ = 1, получаем граничные условия для базовой задачи переходного процесса
Гу = -1; Уу=0, / = 2,3,...,#,-, j = k■, 1
Уу - 0, / = 1, 2, 7 = 1, 2, #г,
(7)
(г0,) =0, 7 = 1,2,...,#,, у**.
еЧ
(7а)
Здесь #2 — число сечений по размаху крыла; Nj — число разбиений в 7-м сечении; 7, к = 1, 2, ...; #2 — номера сечений. Решая задачу с граничными условиями (7), определяем базовые распределения циркуляций диполь-ной системы крыла Гд..
Поясним некоторые особенности решения базовой задачи с граничными условиями (7а). Распределение скоростей перед кромкой и углы атаки сечений крыла связаны между собой линейной зависимостью
=^д(а1]к.(азЪ'’ к = 1> 2, •••> Матрица влияния д(Г0)л/е(а5)у 7=1 ]
не зависит от времени и вычисляется предварительно при последовательном изменении углов атаки сечений на величину .. Обращение мат-
рицы позволяет найти распределение углов атаки сечений по заданным значениям скоростей перед кромкой. Здесь под термином «углы атаки»
д(а$).
понимается распределение (а^. = X Ыу У = 1» 2, ..., #г.
/=1 '
С учетом (7а) получаем базовые значения углов атаки (крутки) сечений эквивалентного крыла. В векторном обозначении они имеют вид
[(“Д]еЧ=^“5)/^°)*,гдей5=(а$)у’*’ •/’ = 1, 2’ •••’ Ы*-
Зная крутку сечений, легко найти базовые распределения циркуляций
вихрей эквивалентного крыла (Г*) . Напомним, что циркуляции вихрей
ея
эквивалентного крыла равны, по определению, циркуляциям присоединенных вихрей рассматриваемого крыла.
Граничные условия (6) и (6а) отличаются от базовых (7) и (7а) множителем У^п~^, поэтому текущее распределение циркуляций и крутку сечений эквивалентного крыла находим как суперпозицию базовых распределений:
к=1 ’
Гч'І^Ь'оІ4. » = 2,3, .... (8а)
к=1
<“*)“ = Д^Ц/о'Т4. "=2>3> - • к—1
7. Скорости перед кромкой крыла определяем в следующей последовательности. На первом временном интервале задаем граничные условия
(1) и находим циркуляцию дипольных систем. Продлеваем систему контрольных точек в область перед крылом и вычисляем начальное поле скоростей гФ, г = О, 1, ..., 7 = 1, 2, , Ы2. Знак минус перед индексом і
означает, что точка находится не на крыле, а в потоке. Далее определяем поля скоростей от базовых дипольных систем У_у/( (г = 0, 1, ..., у,
к = 1, 2, ..., Ы2) или от вихрей эквивалентных крыльев с учетом вклада от
набегающего потока {У-ук) . На этом предварительные расчеты, не со-
еч
держащие временной переменной, количество которых равно Ыг + 1, заканчиваются.
Учитывая вклад начального поля и полей дипольных систем, а также их снос потоком, получаем искомое распределение скоростей перед кромкой крыла:
л-1 N.
иш2' *> - - г-
т=\к=\
Используя (4), можно перейти к полю эквивалентного крыла. Подставляя выражение Ут_п+Ш = {Ут-п+х-Ут_п).к
(7а) получаем
1 &- <-п,1 - (9о)
т=2к=\ 4 *=1 У }
п = 2, 3, ... , 7 = 1, 2, ... ,Ыг,
в уравнение (9), с учетом
ея
гле Лу(т) — У^
гдеАЧ* ~ 0,к У0,к ■
Матрицы в уравнениях (9) и (9а) имеют треугольный вид, их элементы определяются последовательно шаг за шагом по времени, не требуя опера-
ций обращения. После расчета распределения скоростей перед кромкой находим циркуляции вихрей дипольной системы (8) или вихрей эквивалентного крыла (8а),- вычисляем аэродинамические коэффициенты и крутку сечений эквивалентного крыла.
8. Примеры расчета. Выше были описаны два способа построения переходного процесса. Первый предполагал использование дипольных систем крыла, второй — обычную систему П-образных вихрей эквивалентных крыльев. По циркуляциям эквивалентного крыла (8а) можнр сразу определить аэродинамические коэффициенты, в то время как из циркуляций дипольной системы (8) необходимо выделять присоединенную часть вихрей
[2]. Поэтому все практические расчеты, приведенные ниже, проводились по схеме эквивалентных крыльев, т. е. с использованием уравнений (6а) — (9а). ;
На рис. 3 приведен пример расчета прямоугольного крыла с удлинением X = 2,5. Полуразмах крыла делился на восемь равных частей (#, = 8}, хорда — на семь частей (Ых -1). На рис. 3,а показаны относительные
(в долях а) значения местных углов атаки (крутки) сечений эквивалентного крыла в различные моменты времени. Значение I = 1 соответствует прохождению потоком расстояния, равного длине сечения в центральной полосе крыла. Видим, что эквивалентное крыло восстанавливает свою плоскую форму. На рис. 3,6 дано сравнение с результатами работы [2]. Штриховой линией показан уровень стационарного значения су/а. Приведем таблицу
расчетных значений переходной функции для су (отнесенной к ее стационарному значению).
г 0,14 0,71 1,29 1,86 2,43 3,00 3,57 4,14
0,812 0,895 0,936 0,959 0,973 0,981 0,986 0,990
Все полученные выше результаты соответствуют равномерному распределению вихрей по оси X. Если крыло имеет сужение, то обычно применяют разбиение на равные доли хорд. В этом случае геометрическая и временная сетки не совпадают и используется линейная интерполяция для циркуляций [2]. Чтобы не терять преимуществ неравномерных сеток, можно распространить полученные результаты на этот тип сеток.
На рис. 4 представлена информация о треугольном крыле с удлинением Х = 2,5 в тех же разбиениях (#г = 8, = 7). Здесь использовалось разбие-
ние на равные доли хорд. Координаты контрольных точек перед крылом X (в долях хорды корневого сечения), в которых определялось поле скоростей, вычислялись по формуле
Ближайшие к кромкам точки сохраняли топологию контрольных точек крыла, а следующие ряды располагались с равномерным шагом, равным 1/7 хорды сечения с индексом 7 = 1. Никаких поправок на несовпадение временной и геометрической сеток (интерполяций и т. д.) не вводилось. Представляет интерес крутка эквивалентного крыла. Корневая хорда слабо изменяет свой угол атаки во времени в отличие от концевой хорды, которая в начале переходного процесса имеет угол атаки, близкий к нулю. С течением времени эквивалентное крыло становится плоским. Переходная функция изменяется достаточно плавно (см. таблицу).
Г 0,13 0,67 1,21 1,74 2,28 2,81 3,35 3,89
0,804 0,912 0,959 0,979 0,989 0,993 0,996 0,997
Приведенные примеры подтверждают возможность расчета переходных функций по информации о скорости в окрестности передней кромки крыла, при этом размерность задачи построения переходного процесса сокращается. Например, пространственное аэродинамическое нагружение крыла во времени может определяться через крутку эквивалентного крыла.
9. Примеры расчета двумерных течений. В случае плоской пластины вместо (9а) получим следующую систему уравнений:
= Л-I » = 2, 3,
ш=2
(10)
Система уравнений (10) имеет треугольную матрицу с постоянными элементами по диагоналям (когда п -т = const). Элементами матрицы являются компоненты базового поля скоростей (^/)eq. В правые части входят компоненты начального поля РФ.
Для определения базового поля скоростей распределим на пластине систему дискретных вихрей и контрольных точек с безразмерными координатами:
^/=^[V4+0'-1)]; ^=^[3/4+0-0], hj=\,2,...,N.
Как отмечалось выше, начиная со второго временного интервала, эквивалентный профиль становится плоским. Для определения циркуляций вихрей (отнесенных к Vx и хорде профиля) имеем следующую систему уравнений:
N
£|>Л=а; д»°1+2 1-jY '-J='’2.............N-
В работах [4], [5] получены формулы для коэффициентов обратной матрицы этой и других линейных систем, а в работе [6] — выражения для циркуляций:
gN
r/=7yV-r о1)
где
nYl + J-VnYl + J^l, если / = 2, 3, ...,iV-l, m=l V 2mJm=i\ 2т Г ’
nYl + J-\ если i = l,N.
m=l v 2mJ .
(12)
N
Суммируя циркуляции (11), получаем су= 2^Гу = 2л • а.
У=1
С учетом координат вихрей щ = 2^Т^Х-р - = 0, где ХТ = 1/4.
У=1
В работах [6], [7] исследовалось поле скоростей перед кромкой. С учетом вклада набегающего потока в точках Х{ = -(1/4 + /) имеем
Гг/+1
у.=а.2К±2±—I— I - о 1
1 2/ + 1 ^ЛГ+1» 1 1...
Выделим
Уп = а-
2 N
N
(13)
Выполняя условие (7а) (К0)в =1, получаем базовое значение угла атаки
еЯ
к,
N
аеч = -~-тт-. Используя в (11) вместо а величину аеч, находим базовые зна-
ея 2Ы
чения:
-п. г**
(Г,) =— —, / = 1,2,...,#,
V "еЧ N 2# 2/ — 1
(с„) = 2к~, {т2) =0,
N
еЯ
2#
еч
(у.) ,=0 ,
.1 ХТ -N + 1 ’ ’ ’
(14)
,еЯ 2/ + 1 2N К{
Определим начальное поле скоростей. При п = 1 на профиле возникает дипольное распределение вихрей и система уравнений имеет вид
£2>Л=а; ,ч21 '-»-2..
>=1 1 -4(/ - у)
В работе [5] получены выражения для коэффициентов обратной матрицы системы уравнений, соответствующей безциркуляционному обтеканию пластины, в [6] приведены решения для некоторых типов правых частей. В данном случае имеем
N
8 N
Координата центра давления определяется формулой
' 1 ^
Х(1)=_га0)/с<1) + 1/4 = 1
СР 2 / у I 2
1—
2 Ы)
Поле скоростей перед профилем от диполей с учетом составляющей от набегающего потока описывается соотношениями
уО) =п.Е+11+І-В.
-* 2/41 ^N+1+1
7+1
Таким образом, решение при п = 1 получено. Подставляя (14) и (15) в (10), последовательно (для п = 2, 3, ...) определяем величину скорости перед кромкой У^ ■ Условия (6а) и (7а) отличаются множителем У^п~^\ Поэтому текущие значения циркуляции и аэродинамических коэффициентов пропорциональны их базовым значениям:
т^= 0, / = 1,2,...,#, « = 2,3,....
(16)
Отнесем величину У^ к ее стационарному значению (13) и обозначим это отношение как переходную функцию
/^=Г0(и)/к0, « = 1,2,.... (17)
Распределения и отнесем к их значениям при и-»со. В этом случае, согласно (16), имеем
р(П) = у(и-1) . £<!») = /М (18)
Разделив обе части (10) на величину Ко, с учетом (13)—(15), получим систему уравнений для определения переходной функции, дополненной значе-
нием из (17):
,(1) _ N + 1 .
у 2# + 1’
п гп-т+1 р-п
уг*2Ы+ 2п-2т Л1 ^Лт) _ 2# «-1 Л1 ~
^ 2п-2т+\ кМ+п-т^ 2Ы +1 2и-1 к^+п ’ п ^
т=2 1 1
Здесь и далее Д/’^ = Если дополнить последовательность
значений для « = 1, 2, ... величиной -1/2, то получим следующую систему уравнений:
П ТУ-П-М+1 ггП
у2# + 2и-2/и л1_______________________________________________________________________________А Ат) _ 1 А1 , 2 Ппч
^ 2п-2т+\ кМ+п-тч 2 кМ+п’ п а ^
т=1 А1 А1
Уравнения (19) совместно с (12) и (18) представляют решение двумерной задачи в замкнутом виде. Безразмерное значение текущего времени соответствует началу временного интервала
,(„) = «^!5 /7 = 1, 2, .... (20)
Приведем таблицу значений &п\ рассчитанных по формулам (18), (19), в сравнении с точными значениями [1] при числе вихрей N = 4.
І 0,25 0,50 0,75 1,00 2,00 4,00 20,00 50,00
Утеор 0,5556 0,6006 0,6378 0,6693 0,7579 0,8491 0,9702 0,9890
'расн 0,5555 0,6004 0,6376 0,6691 0,7578 0,8490 0,9702 0,9890
Как видно, соответствие с точными результатами оказывается очень хорошим.
10. Вход профиля в ступенчатый порыв. Пусть в момент времени 7 = 0 плоский профиль входит в вертикальный порыв. Зададим безразмерную скорость порыва V- а. На первом временном интервале под воздействие порыва попадает первая контрольная точка. Граничные условия в контрольных точках системы диполей имеют вид
У^+а^О; К(1) = 0; / = 2, 3, N или (^)^=а,
т. е. относятся к типу (6) или (6а). В дальнейшем тип граничных условий сохраняется:
«='. 2,....
где Уц0'* = а.
Следовательно, начиная с п = 1 и далее, на профиле реализуется распределение циркуляций, пропорциональное базовому. Структура системы уравнений для определения скорости перед профилем от дипольных систем имеет вид (9) с начальным полем И9 = а:
В=1-2.......
т=0
(21)
Подставляя, согласно (4), значение Ут-п+\ =(Ут-п+і -Кот_„) в урав-
СЧ
нение (21), переходим к полю эквивалентного профиля:
Ж-,) Д^-а, л-0,1,..., ^'=0.
1=0
(22)
Разделив обе части (22) на Ко, с учетом (13), (14) получим следующую систему уравнений для определения переходной функции (17):
у 2ІУ+2я-2т Щ -------д/т) = 1, и=0, 1, ..., /_1)=0. (23)
^ 2п-2т+\ к*+п~т
Уравнения (23) совместно с (12) и (18) представляют решение двумерной задачи о входе плоской пластины в порыв.
Приведем таблицу расчетных значений в сравнении с точными
данными [1] для числа вихрей N = 5. Дискретные значения времени, в отличие от (20), соответствуют моментам прихода порыва в контрольные точки профиля № ={п-1/4)/Ы.
/ 0,15 0,35 0,55 0,75 1,75 2,75 3,75 4,75
Утеор 0,2406 0,3564 0,4339 0,4928 0,6672 0,7570 0,8120 0,8488
-раем 0,2461 0,3580 0,4346 0,4932 0,6672 0,7569 0,8120 0,8488
Как видно, соответствие с точными данными является вполне удовлетворительным.
11. Связь с интегральными уравнениями. Перейдем к пределу (при N -> оо) в уравнениях (19), (20). Используя формулу Стирлинга
п\— Ппе~п л/27177^1 + у^+- •
можно (см. [4], [6]) выделить главный член в разложении (12) при N>>1:
К? =2$+-.
(24)
Применить это разложение для определения величины К£ от+1 в уравнении (19) при малых значениях параметра (п-пн-1) нельзя. Поэтому выделим часть суммы в (19):
£П—т+1 т=п—п 1
2Ы + 2п-2т к\ 2п-2т + \ ^-л'+я-т
А/
(т)
Зададим значение Я, удовлетворяющее неравенству « п « N. Тогда
±
Я
■л—т+1
^ _ 2п-2т+\
1 т—п—п
М £1
Можно показать, что У 1 = . Обозначив / = п - т +1,
<=1
4^ 2г -1 V N с!х\-----
получим
где
Применяя разложение (24) к остальным слагаемым (19) и переходя к безразмерным значениям времени 1п = п/И, т т=т/Ы, приходим к следующему интегральному уравнению:
Выполняя аналогичные разложения для (23), приходим к интегральному уравнению
Уравнения (25) и (26) имеют одинаковые ядра, но различные правые части. Интегральное уравнение с таким же ядром исследовалось в работах [5], [8]. В качестве искомой функции рассматривалась циркуляция вихрей следа. Интересно отметить, что классические решения задачи о внезапном изменении угла атаки пластины и о ее входе в порыв были получены методами ТФКП и описаны в обзорной работе [1]. Следовательно, решения интегральных уравнений (25) и (26) относительно переходной функции /(О известны. Это — функция Вагнера и функция Кюсснера.
1.Некрасов А. И. Теория крыла в нестационарном потоке. — М.: Изд-воАН СССР. — 1947.
2. Белоцерковский С. М., Скрипач Б. К., Табачников В. Г. Крыло в нестационарном потоке газа. — М.: Наука. — 1971.
3. Г л у ш к о в Н. Н. О сокращении размерности задач нестационарной аэродинамики при малых дозвуковых скоростях // Труды ЦАГИ. — 1996. Вып. 2622.
4. Лифанов И. К., Полонский Я. Е. Обоснование численного метода дискретных вихрей решения сингулярных интегральных уравнений // ПММ. —1975. Т. 39, № 4.
5. Белоцерковский С. М., Лифанов И. К. Численные методы в сингулярных интегральных уравнениях. — М.: Наука. — 1985.
6. Г л у ш к о в Н. Н. Некоторые точные решения дискретной вихревой задачи / В сб.: Проблемы оптимального аэродинамического проектирования // Труды ЦАГИ. — 1977. Вып. 1842.
(25)
(26)
О
ЛИТЕРАТУРА
7. Г л у ш к о в Н. Н. О точности расчета аэродинамических характеристик тонких крыльев и профилей методом дискретных вихрей // Ученые записки ЦАГИ.—1982. Т. XIII, №3.
8. С е д о в Л. И. Плоские задачи гидродинамики и аэродинамики. — М.: Наука. —1980.
Рукопись поступила 10/ІУ1996 г.