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

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

Том XIX 19 8 8

№ 5

УДК 533.6.011.35:533.695

ПОСТРОЕНИЕ ПОВЕРХНОСТИ КРЫЛА В НЕВЯЗКОМ ОКОЛОЗВУКОВОМ ПОТОКЕ ПО ЗАДАННОМУ РАСПРЕДЕЛЕНИЮ ДАВЛЕНИЯ

А. Е. Осовский, В. И. Савицкий

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

В последние годы численные методы невязкого околозвукового течения получили значительное развитие, в результате чего стало возможным их применение к проблеме проектирования формы аэродинамических тел с заданным распределением давления в плоском [1, 2] и трехмерном [3—5] потоке газа.

В работе [3] для решения обратной задачи аэродинамики крыла конечного размаха при околозвуковой скорости набегающего потока используется метод интегральных соотношений, а в работах [4, 5] развивается примененный в [2] подход, состоящий в замене граничных условий Неймана условиями Дирихле.

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

Уравнение изоэнтропийного течения невязкого газа записывается в рамках околозвуковой теории малых возмущений [6]:

-'4'Р хх + В (? у у + ^гг) — срг <р Х2 — 0, (1)

где

А = 1 — - (1 + 7) <р^1 + _1_ <р^ ;

В — 1 + (1 — ч) Мот ’

х, у, г — декартовы координаты (х — в направлении набегающего потока, г — по размаху крыла), — число Маха невозмущенного потока, у — показатель адиабаты, «р (х, у, г) — потенциал возмущения.

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

*/|у=±о = -^т—“• (2)

здесь а — угол атаки крыла, функции г/+(*, г) и у~(х, г) задают координаты соответственно верхней и нижней поверхности крыла. Помимо этого потенциал удовлетворяет условиям непротекания на поверхности фюзеляжа, условиям затухания возмущений на бесконечности и условиям на поверхности вихревой пелены, сходящей с задней кромки

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

еР 1у = ±0 = - 2?х - (I - ЛО - ср2 -<р2. (3)

Для удобства проведения расчетов осуществляется переход к новой неортогональной системе координат 0£г||, позволяющей отобразить бесконечное физическое пространство на конечную расчетную область К I < 1, 1 > ^ > 0, 1 с | < I. Внутри этой области вводится равномерная расчетная сетка и для вычисления поля потенциала возмущения применяется метод верхней последовательной релаксации вдоль линий.

Рассмотрим задачу нахождения функций г/+(х, г) и у~{х, г), для которых выполняется соотношение

ср(х• У)1у = ±о = ср (*■ г>- (4)

где с+ (х, г) и с~(лг, г)—функции, задающие требуемое распределение давления со-ответственно по верхней и нижней поверхностям крыла.

Искомая форма поверхности крыла должна иметь заданную толщину задней кромки

/(*)=/*(*), |

Нг) = У+ К. в, г) — у~(х3 к, г), ) (5)

где х3. к (г)—координата задней кромки, 1*(г)—задаваемое распределение толщины задней кромки.

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

Будем искать распределение потенциала в плоскости крыла, используя формулы (3), (4). Замена этих соотношений их конечно-разностными аналогами приводит к системе квадратичных алгебраических уравнений, разрешая которые последовательно относительно потенциала возмущений получим

-а‘ЦУ (ар)2 - аЦк4к

Чцк =---------------уь-----------’ г=,гп. к> —> гз. к, } —/ф + Ь •"> ■/*> (6)

а\

здесь индекс г'п. к соответствует точкам на передней кромке крыла, 13. к — на задней кромке, /ф — задает граничные узлы, попадающие внутрь фюзеляжа, д — узлы на конце крыла, а индекс к соответствует точке на поверхности крыла. В выражения для входят значения коэффициентов преобразования координат в узлах сетки и известные значения потенциала возмущений в прилежащих узлах, причем при /=/ф + 1 для их определения дополнительно привлекается условие непротекания на поверхности фюзеляжа. Значения потенциала на передней кромке крыла считаются известными, они выбираются из условия заданной толщины задней кромки сечений крыла. В выражения для коэффициентов в формуле (6) входит также квадрат производной Фу в узлах сетки на поверхности крыла. Согласно околозвуковой теории малых возмущений член <Ру в (3) мал по сравнению с (р*, поэтому при расчете значений потенциала возмущений на поверхности крыла первоначально задается величина ц>у, рассчитанная для исходной геометрии крыла, а впоследствии она уточняется.

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

у(х, у=±0, г) = у±(х, г), (7)

где ф+(я, г) и ф~(х, г)—функции, значения которых вычисляются согласно (6). Для решения задачи, как и в предыдущем случае, используется метод верхней последовательной релаксации по линиям. Проведенные численные расчеты показали, что при осуществлении релаксационной процедуры граничные условия Дирихле дают более сильные возмущения, чем граничные условия Неймана. В связи с этим, во избежание появления численной неустойчивости необходимо специально подбирать оптимальную величину коэффициента верхней релаксации. В работе [6] коэффициент верхней релаксации £2 переменный по узлам сетки

., 1 . + ^ . 1 Л V х1+~2] х1~Т])

2 =-------*---------^-------> (8)

ЪхЦ ? !+-£-+!

х!-т1 “

ш

1,0 -

;о”

20°

ж

ЧП°

Рис. 1

где — коэффициент преобразования системы координат, а м — постоянное значение коэффициента релаксации на равномерной расчетной сетке в физическом пространстве.

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

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

ЗуМ _+ ёхкУцк±\ + Ш%кЧцк + .

)Ч ~ 2ВМ

дх

(9)

Л] к

здесь Ц'т определяются коэффициентами преобразования системы координат, взятыми в узлах сетки.

По вычисленному согласно (9) распределению наклона поверхности строится форма профилей в сечениях крыла. Для этого применяется кубическая сплайн-функция, которая на каждом интервале между двумя расчетными узлами задает координаты поверхности по формуле

у(х) =(-&.) (х - х.) + а,

дх и

(х-х& щ (х - х,у

6 Лг

+

+

+ °/+1

(х

■*/)3 + “>1 (X —

6Л;

-Мг. х1+1>х>х1г

(10)

где XI — узлы расчетной сетки, в точках известен наклон поверхности ду ' - х>

(—V

\ дх /г

Лг = х,

<+1

Для нахождения коэффи-

циентов и <1-1 решается система линейных уравнений с тридиагональной матрицей.

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

’к

2 АЧ А,?Ит = 1* (*) —1 (**)• АЧ = ’

1Чф+1

I = ^и. к> 2* — Уф 1, ..., Ук,

где индекс т соответствует точке на поверхности крыла. Коэффициенты матрицы Ац находятся численно по формуле

АЧ = ------>

(И)

(12)

где йу?—малое изменение уут, й/// — соответствующее ей приращение / (гг) после N шагов релаксационной процедуры, /* — приращение при отсутствии возмущений. Итерационный процесс (11) и (12) хорошо сходится уже при Л/=3.

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

Во втором примере рассматривается околозвуковое обтекание комбинации крыло— фюзеляж, типичной для современного пассажирского самолета, при Моо = 0,8 и су = 0,5. С помощью разработанного метода решается смешанная задача аэродинамики: распределение давления на верхней поверхности крыла от борта фюзеляжа до мес-

Ср на построенном крыле.--------------Ср заданное

Рис. 2

п

Верхняя поверхность крыла Нижняя поверхность крыла

Рис. 3

та излома задней кромки модифицировано так, чтобы устранить наличие двух областей сжатия в окрестности х=0,2 и л:—0,6; на этой части крыла определяется форма его поверхности. Модифицированное распределение давления выбирается таким образом, чтобы относительная толщина профилей в сечениях крыла сохранялась неизменной. Для этого, по возможности, поддерживается заданный исходным распределением давления уровень средней линии эпюры (Ср -)- с~)/2, который. согласно приближенной теории тонкого профиля определяет его толщину. Другая часть поверхности от места излома до конца крыла остается неизменной — на ней задается условие непротекания (2).

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

На рис. 5 также представлена форма профилей исходного и построенного крыла. Сравнение приведенных результатов показывает, что желаемого эффекта удалось достигнуть за счет малых изменений формы профилей сечений крыла. На рис. 6 пока- ' заны изменения относительной толщины и крутки, полученные в результате применения предложенной методики решения смешанной задачи обтекания крыла. Распределение толщин профилей сечений крыла по его размаху практически осталось неизменным. Крутка построенного крыла в сечениях 0,1<г<0,4 меньше, чем крутка исходного крыла, при этом величина Ае не превышает 1°. Проверка полученных результатов с помощью программы расчета околозвукового обтекания крыла на основе численного интегрирования полного уравнения для потенциала [10] показала, что коэффициент сопротивления построенного крыла при постоянном значении коэффициента подъемной силы меньше на 8%, чем коэффициент сопротивления исходного крыла. Этот результат объясняется уменьшением волновой составляющей сопротивления в случае расчета обтекания крыла, построенного с помощью изложенного метода.

------исходное крыло;-----------ср модифицированное; --------построенное крыли

Рис. 5

ЛИТЕРАТУРА

1. Шагаев А. А. Построение контура профиля по заданному распределению давлений в сжимаемом потоке газа. — Труды ЦАГИ, 1978, вып. 1925.

2. Т г a n е n Т. L. A rapid computer aided transonic airfoil design method. — AIAA Paper 74-501, 1974.

3. Takanashi S. Iterative three-dimensional transonic wing using integral equations. — Journal of Aircraft, 1985, N 8.

4. H e n n e P. A. Inverse transonic wing design method. —

Journal of Aircraft, 1981, vol. 18, N 2.

5. S h a n k a г V., Malmuth N. D., Cole J. D. Computational

transonic inverse procedure for wing design. — AIAA J., 1982, vol. 20, N 8.

6. Савицкий В. И. Расчет невязкого околозвукового обтекания крыла с фюзеляжем. — Ученые записки ЦАГИ, 1983, т. 14, № 2.

7. Li gh thill М. J. A new method of two-dimensional aerodynamic design. — R&M 2112, April, 1945, ARC.

8. Рыжкова М. В., Серебрийский Я. М. Теоретическое исследование формы тонкого стреловидного крыла с постоянным распределением нагрузки по размаху. — Технические отчеты ЦАГИ, 1953.

9. Седов Л. И. Плоские задачи гидродинамики, аэродинамики.— М.: Наука, 1980.

10. Jameson A., Caughey D. A. Numerical calculation of the transonic flow past a swept wing. — N.-Y. University ERDA Report COO 3077-140, 1977.

Рукопись поступила 14/1V 1986 г.