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

_________УЧЕНЫЕ ЗАПИСКИ ЦАГИ

Том XXIX 199 8

№3-4

УДК 629.735.33.015.3.025.035

МЕТОД РАСЧЕТА ВЯЗКОГО ОТРЫВНОГО ОБТЕКАНИЯ СИСТЕМ КРЫЛОВЫХ ПРОФИЛЕЙ

А. В. Волков, С. В. Ляпунов

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

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

Методика расчета, описанная в настоящей работе, основана на зональном подходе (учете вязкости в тонких пограничных слоях и следах при больших числах Яе и сращивании решений в невязкой и вязкой областях, так называемое вязко-невязкое взаимодействие). Она позволяет определить аэродинамические характеристики систем крыловых профилей (подъемная сила, в том числе максимальная подъемная сила, лобо-

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

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

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

2. Присутствие толстых турбулентных пограничных слоев во взлетной конфигурации и развитых отрывных зон в посадочной конфигурации является неотъемлемой чертой картины обтекания любой механизации.

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

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

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

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

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

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

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

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

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

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

3. Возможные подходы к решению задачи. Можно выделить три подхода к решению задачи об отрывном обтекании системы профилей вязким газом. Наиболее общий и самый сложный подход — это решение осредненных по времени уравнений Навье — Стокса с адекватной моделью турбулентности. В рамках этого подхода все сложные эффекты обтекания могут быть учтены. В настоящее время существует ряд расчетных методов, которые тем или иным образом решают уравнения Навье — Стокса [3], [4]. Однако такие методы требуют больших дорогостоящих ресурсов суперкомпьютеров. Кроме того, несовершенные существующие модели турбулентности не всегда позволяют получать результаты, удовлетворительно согласующиеся с экспериментом.

Один из простых и наиболее старых подходов основан на полуэм-пирической модели отрывного обтекания, когда из предположения о постоянстве давления в отрывной зоне выстраивают форму отрывной области (см., например, [5]). Обычно такие методы могут давать правдоподобные результаты в довольно узком диапазоне чисел Ые и типов взлетно-посадочной механизации. Кроме того, такой подход ориентирован на рассмотрение отрывов в диффузорной части профиля и не приспособлен для рассмотрения обтекания с отрывами иной природы, в частности с отсоединенными отрывами.

Третий подход базируется на решении задачи о вязко-невязком взаимодействии потоков (зональный подход). В настоящее время разработано большое количество численных методов, тем или иным образом решающих эту задачу. При таком подходе предполагается, что все вязкие силы заключены в тонком пристеночном слое и следе. Это позволяет разбить расчетную область на две части: невязкую, где решаются уравнения Эйлера или аппроксимирующие их уравнения, и вязкую, где уравнения Рейнольдса, благодаря тонкости слоя, также значительно упрощаются (уравнения пограничного слоя Прандтля). Решение во внешней (невязкой) и внутренней (вязкой) областях осуществляется последовательно с учетом взаимного воздействия одной области на другую. Иногда утверждается, что такая схема Прандтля учета влияния вязкости позволяет рассматривать лишь безотрывные течения со слабым вязконевязким взаимодействием. Такое представление связано с тем, что уравнения пограничного слоя при заданном внешнем градиенте давления имеют особенность в точках отрыва и присоединения потока (особенность Гольдштейна), не позволяющую интегрировать эти уравнения за точку отрыва. Однако работы [6]—[8] расширили такое ограниченное понимание. В многочисленных работах различных авторов было

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

4. Используемая система уравнений движения вязкой жидкости и газа. Рассматривается течение вязкого сжимаемого газа при больших числах Рейнольдса. При этом воздействие сил вязкости на формирование течения существенно лишь в пределах некоторого тонкого слоя толщин 5, прилегающего к телу, и вязкого следа за ним. Отметим, что течение в вязком пристеночном слое в основном является турбулентным. Запишем в приближении тонкого слоя в криволинейных координатах осредненные по времени уравнения Навье — Стокса для сжимаемого турбулентного потока:

5;(р«) + -^-[0 - М)(Р?)] = о;

95 'г ' дп

8 ( 2) д —[ри + —

95 дп

5-[(1 - л*л)т] -

^-(рн2) + -^-[(1 - Л^лХрио)].- ку,рцо +

др

95

дп1

д ( ч 9 — (рии) + — 95 дп

= _А

дп

(1 - к^п)(ри2) (1 - к„п)ри'2

+ кщри + (1 - ку^п)

лдр

дп

-*т-г(ри'и'\ - к^ри'2\, 95 '

(1)

ди —— т = ц----ри гг.

дп ч

Здесь 5 — длина дуги образующей, п — нормаль, к№ — кривизна поверхности, р — плотность, и, V — тангенциальная и нормальная составляющие скорости, р — давление, р. — коэффициент вязкости, и', V' — пульсационные составляющие скорости, черта над переменными означает осреднение. Граничные условия непротекания на поверхности обтекаемого тела и, и = 0, и стремления возмущенной скорости к нулю при бесконечном удалении от тела формируют задачу, об обтекании тела произвольной формы турбулентным потоком вязкого сжимаемого газа.

Поток, описываемый системой уравнений (1),- в дальнейшем называется реальным вязким потоком (РВП). В работах [6], [7] предложено моделировать свойства РВП так называемым эквивалентным невязким потоком (ЭНП), описываемым системой уравнений (1) с нулевой правой частью (уравнения Эйлера).

ю

Граничные условия на бесконечности аналогичны граничным условиям РВП. Однако граничные условия на поверхности обтекаемого тела подбираются из условия совпадения свойств РВП и ЭНП на расстояниях больших, чем 8 от обтекаемого тела. Запишем разницу системы уравнений ЭНП (переменные этой системы уравнений имеют индекс /') и системы РВП. Произведя интегрирование по нормали от 0 до 2

5, разделив на р^и,*, и учтя эквивалентность потоков при я >8, получим [7] (здесь символ означает значение величины на стенке):

Р™

1 а

Р *

(р/ч^е) + —^ ~ + (н + 2 -

ит (Ьа$ и т аз

— с/ + ку,($п + сТ) +

2

Ар

(л - рМп + |Ри'2<1п

____________ 1 с1

Р/и>ы/и' Рти'т ^

-С/+Е{'

(р- — + К (е + 8*) + —| (р1 - р)йп

РпЛ'уу

РЫ/ит о

}р и'

2, 1 (1 ап +

О

| (-р и'ю')йп = + 8*) + Е2-

Здесь М — число Маха на стенке в ЭНП, т. е. М = М,-*,, и введены следующие обозначения:

* 1 5 1 5 8* =--------[ (р,иг- - ри)(1п; 0 =----=- | [р и(иш - и)- р/м,(м/н, - и)]Лг;

Рм,ип» о

•(-ри'о')

АР = Рт - Р*; 0л = —Ц"}(р,-И/Ч - рму)^и; ст = }■ 2

Р/и>ит о 0

СЫ\

С/=Т

(2)

- Р;И1М/Н’

С учетом следующих полуэмпирических соотношений, полученных в работе [9],

§ 5

0„ = (о + 8*)^-: ](л - р)с!п = уР/м,м^*(0 + 8*) +|ри,2(1п-,

о

| ри’2 с1п = 0,05рги,иги, о

5

| ри'у'</я = 0,03рги,«(

Я-.1

Я

о

0; \ ри'2 с!п = 0,12РмЧ2К\ 1

); }р«'2

о

я

0;

Я

</28*

0

И

можно показать в [7], что значения величин Е\, Еі обычно малы. Большая часть входящих в них слагаемых имеет второй порядок малости относительно величин 5* и 0, а те, которые имеют первый порядок малости относительно 8, умножаются на малые коэффициенты. В связи с этим величины и £2 в дальнейшем будут опущены.

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

Нетрудно показать, что и,(я) = ие + 0(8), где ие — значение тангенциальной скорости на внешней границе пограничного слоя (я = 8),

и, следовательно, переменные (2) могут быть приближенно записаны в виде

что делает более понятным их физический смысл: величина 8* аналогична толщине вытеснения пограничного слоя, величина 0 аналогична толщине потери импульса пограничного слоя.

Теперь систему уравнений (3) можно интерпретировать для ЭНП следующим образом. Уравнение (За) представляет собой граничное условие, связывающее величину нормальной скорости на поверхности обтекаемого тела и толщину вытеснения пограничного слоя. Иначе говоря, на поверхности контура необходимо расположить источники, интенсивность которых такова, что линии тока оттесняются на расстояние 8*. Такой способ учета воздействия пограничного слоя на формирование течения был предложен Лайтхиллом [10] еще в 1958 г. и был назван методом эквивалентных источников. В этой же работе показана аналогичность такого способа учета вязкости способу, предложенному Прандт-лем, когда из соображений сохранения расхода в приграничном слое форма тела изменяется на толщину вытеснения пограничного слоя.

Уравнение (36) аналогично уравнению Кармана для пограничного слоя. Это уравнение отражает закон сохранения импульса в пограничном слое в продольном направлении. После определения решения в эквивалентном невязком потоке использование простого соотношения (Зв) позволяет получить решение в реальном вязком потоке.

задачи, описываемой системой уравнений Эйлера со следующими условиями на границе обтекаемого тела:

(3)

Количество неизвестных в системе уравнений (3) превышает количество уравнений. В соответствии с [11] дополним эту систему уравнений следующими уравнениями:

се

(реме//]б);

(4)

0 (1иР

а{

г1/2 1с1!2

СхЕО0 кст

(ит/ие)(Ь/8) Я + Я,

0 (1ие ие йЬ

ЕО

ие с1ь

1 +0,075М2 ( 0 \ 1 + 0,2М

Ч1 + 0,1М2 J

(5)

е _ о*

где Н\ = —-—. Первое из этих уравнений представляет собой опреде-0

ление коэффициента эжекции С£, характеризующего количество газа, эжектируемого вязким течением из внешнего невязкого течения. Уравнение (5), называемое уравнением «запаздывающей эжекции», получено в [11] из уравнения для турбулентной кинетической энергии и представляет, по существу, интегральный аналог дифференциальной модели турбулентности одного уравнения Брэдшоу [12]. В этом уравнении ст —

коэффициент касательного напряжения, ст = ттах/ре«2, касательное напряжение оценено при у/5 = 0,2, с' = йсх/йс£ . Знак Е(2 означает, что соответствующая величина получена при равновесных условиях при неизменных формпараметре и коэффициенте эжекции, ах — параметр структуры турбулентности, ЦЬ — обезразмеренная длина диссипации, ит!ие ~ отношение продольной скорости в точке максимального напряжения сдвига к скорости на внешней границе пограничного слоя.

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

iR.ee, вычисленное по толщине потери импульса.

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

и, вообще говоря, простирающаяся до бесконечности. Вязкий след, как и в работе [11], является продолжением развития пограничного слоя на профиле за заднюю кромку вдоль линии разрыва. В соответствии с этим след разделен на две части: верхнюю (индекс и) — продолжение пограничного слоя с верхней части контура и нижнюю (индекс Г) — продолжение пограничного слоя с нижней части контура. Для моделирования

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

Рш ® Р/7

Здесь и ниже символ «с» означает значение величины в следе.

Кривизна вязкого следа приводит к перепаду тангенциальной составляющей скорости вдоль линии разрыва. В соответствии с соотношением (Зв) разница между значениями скорости в ЭНП и РВП на линии разрыва соответственно на верхней и нижней поверхностях составляет

к„ — кривизна линии разрыва. ;

Разница скорости между верхней и нижней сторонами разрыва составляет

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

5. Методы вязко-невязкого взаимодействия. Зональный подход требует согласованного решения системы нелинейных уравнений для невязкой и вязкой областей потока. Существует множество способов учета вязко-невязкого взаимодействия [8]. Рассмотрим некоторые из них более подробно.

Перепишем уравнения (4) и (За) в виде

Здесь и далее индекс характеризующий параметры ЭНП на стенке, опущен. Так как замыкающие соотношения пограничного слоя выведены для несжимаемой жидкости, вводится формпараметр эквивалентного несжимаемого течения Н, который связан с Я соотношением [11]:

(Ди)с = иш - ии * ку,ис(ви + 0/ + 8* + 8*); ис = ±(иш +,/). (7)

(8)

(9)

Учитывая уравнение импульсов (Зв) и принимая во внимание, что

йНл

сЬ

<т_

Н{ =

СІН

уравнения (8) и (9) можно переписать:

сЕи = -(Н + \)Н1в^ +

ІПс\(ІН . 1 гг

1 ~3я~ 2" ^

и;

(10)

-(Я + 1)(Я-О,4М2)0^

и.

(11)

В «прямом» методе расчета решение внешней невязкой задачи с известным значением величины нормальной составляющей скорости на поверхности V дает значение величины градиента скорости на внешней границе пограничного слоя с1и/с18. Следует отметить, что решение внутренней вязкой задачи может быть получено на основе дифференциального (конечно-разностного) или интегрального подхода. В настоящей работе анализ существующих методов вязко-невязкого взаимодействия производится на основе интегрального подхода, что, вообще говоря, не ограничивает общности рассмотрения. Решение получается после интегрирования системы дифференциальных уравнений — уравнения импульсов (36), уравнения для «запаздывающей» эжекции (5), а также уравнения (10), записанного в виде

0

ст

<1$

Щ

с£-Я,{1с/-(Я + 1)£^}

(12)

Новое значение составляющей скорости V на обтекаемом контуре получается из уравнения (11), которое, после подстановки в него (12), можно переписать в виде

Т & Т

(13)

где

Я

НН{ - Я

Р = -(Я + і)(і-О,4М2г)0; Є=ІС/+£|^; щ=-К 1 <

ч

Область применения прямого метода ограничена лишь безотрывным обтеканием, и при наличии отрыва решение становится невозможным в связи с сингулярностью (известной как особенность Гольдштейна) записанных таким способом уравнений в точке отрыва (значение производной Н{ в уравнении (12) обращается в ноль (рис. 1)). Фурье-анализ устойчивости «прямой» схемы вязко-невязкого взаимодействия, приведенный в работе [13], показал, что в присоединенной области по-

тока вязко-невязкие итерации сходятся к решению при условии нижней релаксации поправок и и расходятся в отрывной области.

Однако уравнения (10), (11) можно переписать и в другой форме, выразив градиенты формпараметра и скорости через значения коэффициентов эжекции се и нормальной скорости и:

Это обратная форма записи уравнений пограничного слоя, так как в этом случае на внешней границе пограничного слоя задается величина V (являющаяся, по существу, производной от толщины вытеснения 5*), а рассчитывается скорость и на внешней границе. Следует отметить, что при этом в процессе расчета сингулярностей не возникает. Впервые это было продемонстрировано в работе [14] для ламинарного пограничного слоя. На основе этого были созданы обратные вязко-невязкие методы расчета обтекания (например, [15]), когда и задача о внешнем невязком обтекании, и внутренняя задача решаются в обратной форме. Задание скорости на границе обтекаемого контура позволяет получить новую форму контура, что, в свою очередь, позволяет произвести расчет пограничного Слоя в обратной форме. Такой метод, однако, не получил широкого распространения, во-первых, из-за трудностей создания обратно методов расчета невязкого обтекания, во-вторых, из-за расходи-

О = ~(НН1 - Я,)(1 - 0,4Мгг)г; % =

и

мости итераций в присоединенной части пограничного слоя, да и в отрывной части скорость сходимости сильно зависела от размеров расчетной области [13].

на профиле

в следе

Применение «полуобратной» схемы вязко-невязкого взаимодействия позволило получить сходящееся решение как в отрывной, так и в безотрывной области течения на основе итерационного взаимодействия прямой задачи расчета невязкого обтекания и

2

н обратной вязкой задачи расчета погра-

Рис. 1. Профиль скоростей и замыкающие соотношения Я, (я),

для турбулентной части пограничного слоя

Присоеди-

ненный

слой

о

2 | 4 6 8 10 12 14

Я

4 6 8 10 12 14

'у! ■ ' ■ '■ ■ і............

Оторвавшийся

спой

£1е..=Ю

ничного слоя. Первые исследования таких схем были выполнены в работах [16]—[18]. Расчет во внешней области течения дает значения величины «невязкого» градиента скорости йи1 !<к на обтекаемом контуре, расчет обратным методом вязкой области дает величину «вязкого» градиента Разность между значениями этих градиентов (или между значениями скоростей) используется для оценки новой величины нормальной скорости V.

Оптимальное значение для коэффициента релаксации со вычислено в работе [13]. Метод сходится как с присоединенным, так и с отрывным пограничным слоем. Однако опыт применения этого метода авторами настоящей статьи показал, что такой способ позволяет получать лишь непротяженные отрывные зоны. Для эффективного корректирования толщины вытеснения требуется гладкость получаемых скоростей, а следовательно, и высокий уровень сходимости в невязком течении. В случае наличия протяженных отрывов и больших толщин вытеснения достижение высокого уровня сходимости влечет за собой значительное увеличение времени расчета. Кроме того, как показано в [13], при достаточно «толстом» относительно шага разностной сетки пограничном слое возникает так называемая «высокочастотная» неустойчивость.

Использование одновременных методов позволило получать решения с протяженными отрывными зонами, с отрывами из-под скачка уплотнения. Суть метода состоит в одновременном, согласованном решении вязких и невязких уравнений, причем во внутренней, вязкой области уравнения записываются в форме, не имеющей особенностей в точках отрыва и присоединения потока. В последнее время было развито большое количество одновременных методов сращивания вязконевязкого взаимодействия [19], [20]. Одновременный метод сращивания на основе метода Ньютона для расчета обтекания профиля вязким сжимаемым газом (код 18Е5) изложен в [21]. В указанных работах демонстрируются решения с локальными ламинарными отрывами и турбулентными отрывами из-под скачка уплотнения. Однако при расчете этим методом приходится обращать довольно большую и полную матрицу, что требует больших вычислительных затрат.

«Квазиодновременные» методы более экономичны. Можно выделить два типа «квазиодновременных» методов. В первом типе уравнение, тем или иным образом аппроксимирующее влияние невязкого течения на характеристики пограничного слоя, решается совместно с вязкими уравнениями. Это аппроксимирующее уравнение выражает связь между толщиной вытеснения и скоростью внешнего потока через интеграл линейной теории потенциальных течений жидкости (например, интеграл Коши). Этот метод был предложен в [22]. Особенности его использования исследовались в [23], [24]. В работе [25] реализовали этот метод в программе расчета трансзвукового обтекания профиля У18ТКАР5. В последней работе авторы демонстрируют расчеты с локальными отрывами из-под скачков уплотнения. Однако при этом не демонстрируются решения с развитыми отрывными зонами в окрестности максимальной подъемной силы. «Квазиодновременный» метод вязко-невязкого взаимодействия по способу [22] был применен также к задаче расчета обтекания многозвенного профиля в работе [1].

В другом типе «квазиодновременных» методов, наоборот, какая-либо аппроксимация влияния пограничного слоя на внешнее течение добавляется к уравнениям внешнего невязкого обтекания в качестве граничных условий. Такой метод был впервые описан и реализован в [26], [27] в программах расчета отрывного обтекания профиля несжимаемым газом. Авторами данной статьи «квазиодновременный» метод был применен в задаче вязкого трансзвукового обтекания профиля [28].

В настоящей работе реализован смешанный метод учета вязконевязкого взаимодействия. В областях ламинарного и присоединенного тонкого турбулентного (Я < 1,4) пограничного слоя используется прямой метод (расчет пограничного слоя по известному распределению давления). На остальной части обтекаемого контура применяется «квазиодновременный» метод учета вязко-невязкого взаимодействия. Решение внутренней вязкой задачи производится при заданном значении величины нормальной скорости V. Методом Рунге — Кутта — Мер-сона пятого порядка точности с автоматическим выбором шага производится интегрирование системы дифференциальных уравнений:

1) уравнения импульсов (36);

2) уравнения (14) (являющегося следствием уравнений (За) и (4));

3) уравнения «запаздывающей» эжекции Грина [11] (5).

Рассчитанные параметры пограничного слоя позволяют получить

значения переменных Т, Р, 0 уравнения (13) (также являющегося следствием уравнений (За) и (5)). Уравнение (13), записанное в виде

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

которое следует из уравнения (6), здесь Тс=ТиТі\ Рс = РиТ{ + Р{Ги;

Скачок тангенциальной составляющей скорости (уравнение (7)) в настоящей работе не учитывался.

6. Замыкание системы уравнений турбулентного пограничного слоя.

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

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

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

П-Р^-Ои = 0

ес =йиТ1+д1Ти.

1, 0 < г| < г|*,

^. = 1 + с11п(л)-с2/’(т,); л = у/ь; /?(л)= /(^

где

/)] =-5,25.

+ Д эе

Отметим, что данное семейство профилей скорости включает в себя профили с зоной возвратных течений, характерные для отрывного пограничного слоя (см. профили скорости на рис. 1). Здесь Сг, по существу, определяет скорость возвратного течения в зоне отрыва. Функция 7(5) определяет следовую часть профиля скорости (различные способы задания этой функции рассматривались в работах [6], [8], [26]). Величина л* ограничивает зону возвратного течения по высоте пограничного слоя. В присоединенном потоке она равна нулю, а в отрывном определяется эмпирически через величину а = 8*/д. Как показано в [6], при больших значениях Я (в зоне развитого турбулентного отрыва), (1г\*/<1а «2,3. В промежуточной области, между безотрывным обтеканием и областью с сильным отрывом, связь между л* и а определяется квадратичной функцией:

0; 0 < а < ао,

2,3 / \2

-------\Ча1 -а0) , а0<а<аь

21а] -а0)

2,3а + 1 - 2,3; а] < а < 1.

■ Величина ао ограничивает область отрывного обтекания и определяется выбором функции /(£). Значение величины а однозначно определяет формпараметры:

тг _ а .

~0/8’

Н,-Н (1-11,

(16)

а следовательно, параметрически и зависимость Н\ (я) в (15).

После многочисленных параметрических исследований зависимость Н\{и) была представлена в виде

на профиле: Н\ =

0,63 Я +1

Я+0,75; Я <2,61,

Я -1

, 25,07 35,78.

7,27 + ^==— ----------==—; Я >2,61;

Я-0,42 Я

в следе:

Я,

0,40Я +1 77 ■п ;

4,00

Я -1 1,83

Я <2,61,

5,12.

Я -1,18 Я

; Я >2,61.

(17)

Коэффициент трения су в соотношении (15) выражался в следующем виде [29]:

Остальные замыкающие соотношения в (15) для присоединенного пограничного слоя совпадают с использованными в [11], а в области развитого отрыва (при Н > 5) они равны:

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

В оторвавшемся потоке необходимо также изменить величины аь Ь/8 и мт/ме, которые в присоединенном потоке имели соответственно значения 0,15, 0,08 и 0,67, а в отрывном течении эти величины равнялись 0,16, 0,08(1-л*) и 1-0,7С2 [30].

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

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

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

Описанная ниже методика расчета положения ламинарнотурбулентного перехода была разработана В. М. Галкиным.

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

с/ =

сееоо =0,024С2; сх - 4,4с|.

(18)

выражения: = а = С| + 0,5С2(1 + л*)

ристиками пограничного слоя. Эта связь устанавливается с помощью решения уравнения Орра — Зоммерфельда, в которое подставляются профили скорости из семейства Фолкнера — Скэн с различными значениями формпараметра. На основании этих решений устанавливаются зависимости коэффициентов нарастания для различных частот неустойчивых возмущений от формпараметра и числа Яе. Огибающие интегралов коэффициентов нарастания возмущений аппроксимируются соотношениями, приведенными в [31], требующими решения дополнительного обыкновенного дифференциального уравнения относительно коэффициента нарастания возмущения п.

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

Система интегральных уравнений ламинарного пограничного слоя дает возможность продолжать расчет за точкой ламинарного отрыва. В настоящей программе при продолжении расчета за точку ламинарного отрыва принимается Cf = 0, и расчет устойчивости ведется так же, как и для присоединенного слоя, при этом быстрое нарастание возмущений обеспечивается за счет быстрого роста формпараметра Н. При выполнении критерия перехода за точкой ламинарного отрыва считается, что происходит присоединение турбулентного пограничного слоя. Типичный отрывной пузырь имеет большой градиент скорости на внешней границе ие и толщины потери импульса 9 вблизи точки присоединения. Эти величины на малом расстоянии претерпевают заметные скачки Аие и Д0, особенно при небольших числах Яе. Соотношение, связывающее скачки Аие и Д0, легко получить, если проинтегрировать уравнение импульса по х на малом интервале Ах без учета трения на стенке с у. Это

соотношение имеет вид 4г = -(2 + Н - м2)-—и применяется для оп-

6 ’ ие

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

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

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

При решении задачи о невязком обтекании системы профилей интенсивности источников и вихрей находятся из условия непротекания, выполняемого в центрах каждой панели. После получения решения задачи о невязком обтекании системы профилей за каждым из профилей выстраивается линия следа, соответствующая форме линий тока, сходящих с задних кромок элементов. Правильный расчет параметров пограничного слоя в следе во многом определяет точность решения. След от передних элементов, проходя над нижерасположенными элементами, существенным образом может влиять на формирование течения. В расчетах фиксировалась общая длина следа Лю за каждым из элементов. Расстояние от задней кромки до первого узла в следе, как правило, принималось равным значению длины последней панели профиля. Далее расстояние между узлами монотонно увеличивалось, но не превышало 2,5% хорды основного профиля системы. Общее количество точек в следе N варьировалось от 100 при Лю = 4 до 200 при Луу = 8 хорд.

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

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

ь Для обеспечения заданного перепада нормальной скоро-

сти в следе на панели следовой линии располагаются источники с интенсивностью а*, = — -т-|м8* г.

№ и си 1

В зонах с толстым и оторвавшимся пограничным слоем Н> 1,4 использована «квазиодновременная» методика учета вязко-невязкого взаимодействия, описанная выше. Интенсивности особенностей на профиле и в следе находятся при выполнении следующих граничных условий:

Ть = Р^ + 0ц\ = Р„4±. + (2„и, где Т, Р, 0, Ту,, Рк, (2„ — пара-

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

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

на профиле:

ди

ЙУ

-кю

8у

дп

в следе:

ди

дя

о -

Зи да

дп и дп 1.

здесь индекс 0 означает значение скорости в центре панели за вычетом вклада самой панели. ••

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

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

9. Сравнение результатов расчета с экспериментальными данными. По разработанной на основе описанной методики программе М и ЬТ Г/1 Э были

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

На рис. 2 демонстрируются результаты профиля СА(\У)-1 с 30%-ным закрылком типа Фаулера 53 = 35°. Величина щели вырьировалась от 2 до 5% при постоянной величине перекрытия 1,2% [32]. Испытания проводились при числе Яе = 2,2 млн., М = 0,13 в малоскоростной 213 х 305 см трубе со специальной рабочей частью 213 х 91,4 см прямоугольного сечения. Измерялись как рас-

Рис. 2. Сравнение результатов расчета с экспериментом для профиля ОА(\\0-1 с 30%-ным закрылком, 83 = 40°, Яе 2,2 млн.:

расчет МШЛТУШ эксперимент [32]

------- щель 3% Д

• • • щель 2% +

------- щель 5% О

пределение давления, так и весовые характеристики — подъемная сила, сила сопротивления, продольный момент. На рисунке приводится распределение давления на основном элементе и закрылке при а = 4°. Видно, что в расчете, как и в эксперименте, при увеличении ширины щели с 3 до 5% существенно изменился характер обтекания закрылка от безотрывного (3%) до существенно отрывного (5%). Здесь же приведены зависимости величины подъемной силы по углам атаки для конфигураций с шириной щели 2, 3 и 5%. В левой части рисунка результаты расчетов, в правой — экспериментальные данные. Наблюдается хорошее качественное согласование расчетных данных с результатами испытаний. Некоторое превышение расчетной величины су можно объяснить неучетом пограничного слоя в подрезке основного элемента.

Аналогичные испытания были проведены для этого профиля и с 29%-ным закрылком для 83 = 10, 20, 30 и 40* [33]. На рис. 3 демонстрируется хорошее согласование расчетных и экспериментальных значений для углов отклонения 53 = 10, 20, 30°. Для 83 =40° величина критического угла атаки на градус меньше, чем в эксперименте, и уже при а = 11° наблюдается значительное уменьшение несущих свойств исследуемой конфигурации. Темп приращения величины максимальной подъемной силы в зависимости от угла отклонения закрылка в расчете несколько ниже, чем в эксперименте, тем не менее основная зависимость качественно отражена.

Рис. 3. Сравнение результатов расчета с экспериментом для профиля ОА(\\0-1 с 29%-ным закрылком, 53 = 10, 20, 30, 40°, Яе = 2,2 млн.:

---- расчет МииПУК; 0, Д, V, О — эксперимент [33]

Для угла а = 10° (83 = 40°), являющегося критическим углом атаки, результаты расчета в деталях представлены на рис. 4. Здесь приводятся сравнение расчетной эпюры давления с экспериментом, величина поверхностного трения на основном элементе и закрылке, форма линий тока, а также профили скорости в трех сечениях. В целом эпюра давления удовлетворительно согласуется с экспериментом. Некоторое рассогласование в диффузорной части закрылка можно объяснить большой

Рис. 4. Расчет обтекания профиля ОА(\\0-1 с 29%-ным закрылком,

М = 0,21, Яе = 2,2 млн.; линии тока при а = 10° с зоной возвратного

течения:

---распределение давления в сравнении с экспериментом [33]; О — профили скорости в сечениях Л, В, С

относительной толщиной задней кромки, что в настоящем расчете не принималось во внимание. Расчет, как видно из графика поверхностного трения, показывает наличие непротяженных зон диффузорного отрыва (область отрыва занимает 14% основного элемента и 3% закрылка), что также согласуется с экспериментальными наблюдениями [33]. Несмотря на малую протяженность зон отрыва на основном элементе и закрылке, режим на данном угле атаки является режимом с максимальной подъемной силой. Ограничение в величине максимальной подъемной силы можно объяснить присутствием значительной отрывной зоны в следе основного элемента (на рисунке область возвратного течения обозначена точками). Линии тока, полученные в расчете, огибают протяженную зону отрыва, начинающуюся на основном элементе и развивающуюся над закрылком. Как видно из рисунка, зона отрыва отделена от закрылка слоем безотрывного течения. На этом же рисунке представлены профили скорости в сечениях: А — область задней кромки основного элемента, В — сечение в следе с наибольшей областью возвратного течения, С — сечение в следе далеко от контура. Таким образом, видно, что расчеты по вышеописанной методике моделируют течения с отсоединенными отрывными зонами, наблюдаемыми в экспериментах [34].

Пример расчета интегральных аэродинамических характеристик для профиля с предкрылком приведен на рис. 5 для М = 0,201, Ле = 2,83 • 106, 8пр = 43,9°. Здесь приведены зависимость величины подъемной силы по углам атаки, величина сопротивления и момент тангажа (относительно четверти хорды основного элемента) в зависимо-

Рис. 5. Сравнение результатов расчета с экспериментом для профиля ЫА5А с предкрылком:

----- расчет МииПУК; О — эксперимент [35]

Л

Рис. 6. Сравнение результатов расчета с экспериментом для четырехэлементной системы профилей, образованных профилем NASA:

---- расчет МЧЛЛТУК; О — эксперимент [35]

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

На рис. 6, 7 демонстрируются возможности метода расчета аэродинамических характеристик четырехэлементной системы профилей. На рис. 6 приводятся интегральные аэродинамические характеристики для изолированного профиля, а также для двух конфигураций, полученных из этого исходного профиля. На этом же рисунке приводятся экспериментальные результаты работы [35]. Конфигурация с наибольшими несущими свойствами имела следующие параметры: хорда пред-

крылка Ипр =16%, угол его

отклонения

42, Г

щель между предкрылком и основным элементом

£гпр = 1,82%, хорда дефлектора Лд — 15%, 5Д = 30°, <?д= 1,1%, величина перекрытия Од = -0,32%, хОрда закрылка А3 = 29%,

<fc

\

JrV

Рис. 7. Сравнение расчетной эпюры давления с экспериментом для четырехэлементной системы профилей, образованных профилем NASA, при а = 8°, Яе = 2,83 млн.:

---- расчет МиЬТ1У18; О — эксперимент [35]

83 = 49,7°, (?3 = 1,69%, величина зазора между дефлектором и закрылком составляла 03 = 0,07%. Другая конфигурация имела менее критические параметры: 8пр = 50°, (7пр = 3,56%, 5Д = 0, Ga = 1,94%, Од = -0,01%, 83 = 16°, G3 = 2,44%, 0з = 1,11%. На рис. 7 для первой конфигурации расчетная эпюра давлений на всех четырех элементах сравнивается с экспериментальной.

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

ЛИТЕРАТУРА

; 1. King D. A., Williams В. R. Developments in computational

methods for high-lift aerodynamics // Aeronautical Journal, August, September. — 1988.

2. В r u n e G. W., McMasters J. H. Computational aerodynamics applied to high-lift systems // In: Applied Computational Aerodynamics. — 1989.

3. К о g a n A., M a г о m S. Navier — Stokes calculation of the flow about wing-flap combinations // ICAS-92-4.1.4. — 1992.

4. С а к о в и ч В. С., Сорокин А. М. Применение неструктурированных сеток для расчета вязкого обтекания многоэлементных профилей // Ж. вычисл. матем. и матем. физ. — 1997, № 7.

5. С т е п а н о в Ю. Г., Р ы ж о н к о в М. В. Приближенный метод расчета отрывного обтекания профиля с отклоненным щелевым закрылком // Труды ЦАГИ. — 1983. Вып. 2213.

6. L е В а 11 е u г J. С. Strong matching method for computing transonic viscous flows including wakes and separations // La Recherche Aerospatiale. —

1981, N 3.

7. L о с k R. С., F i r m i n М. С. P. Survey Of techniques for estimating viscous effects in external aerodynamics // Proceedings of IMA

Conference on Numerical Methods in Aeronautical Fluid Dynamics, 30 March — 1 April 1981 / Edited by P. Roe, Academic Press. — 1983.

8. Lock R. C., Williams B. R. Viscous-inviscid interactions in external aerodynamics // Progress in Aerospace Sciences. — 1987. V. 24.

9. E a s t L. F. A representation of second-order boundary layer effects in the momentum integral equation and in viscous-inviscid interactions // RAE TR 81002. - 1981.

10. L i g h t h i 11 M. J. On displacement thickness // J. Fluid Mech. — 1958. V. 4.

11. Green J. E., Weeks D. J., В room an J. W. F. Prediction of turbulent boundary layers and wakes in compressible flow by a lag-entrainment method // RAE Technical Report 72231 (ARC-RM 3791). — 1977.

12. Б p э д ш о у П. Турбулентность. — М., 1980.

13. W i g t о n L. В., H о 11 М. Viscous-inviscid in transonic flow // AIAA Paper. N 81-1003. - 1981.

14. С a t h e r a 1 1 D. and Mangier K. W. The integration of the two-dimensional laminar boundary layer equations past the point of vanishing skin friction // JFM. - 1966. V. 26, p. 1.

15. С a 1 v e r t W. J. ASME Paper 82-GT-135. - 1982.

16. L e В a 1 1 e u r J. C. Viscous-inviscid flow matching: numerical

method and applications to two-dimensional transonic and supersonic flows // La Recherche Aerospatiale. — 1978. " 1

17. С a r t e r J. E. A new boundary layer inviscid iteration technique

for separated flow // AIAA Paper-1450. — 1979.

18. Б p у т я н М. А. Исследование отрывного обтекания симметричного профиля в несжимаемой жидкости // Ученые записки ЦАГИ. — 1978. Т. IX, № 6.

19. G i 1 m е г В. R. and В ;r i s t о w D. R. Simultaneous perturbations to viscous and inviscid equations // AIAA Paper-1239. — 1981.

20. H о u w n i k R. and V e 1 d m a n A. E. P. Steady and unsteady separated flow computations for transonic airfoils // AIAA-1618. — 1984.

21. G i 1 e s М. B., D r e 1 a M. Two-dimensional transonic aerodynamic design method // AIAA J. — 1987. V. 25.

22. Veldman A; E, P. New, quasi-simultanepus method to calculate interacting boundary layers // AIAA J. — 1981. V. 19.

23. Edwards D. E., С a r t e r J. E. A quasi-simultaneous finite

difference approach for strongly interacting flow // Numerical and Physical

Aspects of Aerodynamic Flows III. — New York, 1986.

24. С e b e с i Т., Jau J., V i t i e 11 о D., Chang К. C: Prediction of post-stall flows on airfoils // Numerical and Physical Aspects of

Aerodynamic Flows IV. — 1990.

25. Veldman A. E. P., Lindhout J. P. F., Boer E.,

Somers М. A. M. VISTRAFS: A simulation method for strongly interacting viscous transonic flow // Numerical and Physical Aspects of Aerodynamic Flows IV. - 1990.

26. С г о s s A. G. T. Boundary layer calculation and viscous-inviscid coupling // ICAS-2.4.1. — 1986.

27. Величко С. А. Решение задачи обтекания профиля при

больших числах Рейнольдса // Изв. АН СССР, МЖГ. — 1990, № 4.

28. Volkov А. V., Lyapunov S. V. Numerical prediction of transonic viscous separated flow past an airfoil // Theoretical and Computational Fluid Dynamics. — 1994, N 6.

29. S w a f f о r d T. W. Analytical approximation of two-dimensional

separated turbulent boundary layer velocity profiles // AIAA J. — 1983. V. 21, N6. ••

30. M e 1 n i k R. E. and Brook J. W. The computation of viscous-inviscid interaction on airfoils with separated flow // Numerical and Physical Aspects of Aerodynamic Flows III. — New York, 1986.

31. D r e 1 a М., Giles М. B. Viscous-inviscid analysis of transonic and low Reynolds number airfoils // AIAA J. — 1981. V. 25, N 10.

32. W e n t z W. H., Jr., О s t о w a r i C. Additional flow field studies of the GA(W)-1 airfoil with 30-percent chord fowler flap including slot-gap variations and cove shape modifications // NASA CR-3687. — 1983.

33. W e n t z W. H., Jr., S e e t h a r a m H. C. Development of a fowler flap system for a high performance general aviation airfoil // NASA CR-2443. - 1974.

34. Петров А. В. Некоторые типы отрывного обтекания разрезных крыльев // Ученые записки ЦАГИ. — 1977. Т. VIII, № 2.

35. О ш а г Е., Z i е г t е п Т., Hahn М., S z р i г о Е. Twodimensional wind-tunnel tests of a NASA supercritical airfoil with various high-lift systems. V. II//NASA CR-2215. - 1973.

Рукопись поступила 25/HI 1997 г.