УДК 629.7.072
Л. Л. Теперин1'2, Л. Н. Теперина2, Мьё Тхейн1
1 Московский физико-технический институт (национальный исследовательский университет)
2 Центральный аэрогидродинамический институт им. профессора Н. Е. Жуковского
Проектирование аэродинамической поверхности крыла в компоновке с горизонтальным оперением на трансзвуковом крейсерском режиме полета
На примере компоновки крыла и горизонтального оперения пассажирского самолета излагается методика проектирования формы сечений крыла, обеспечивающих минимальные индуктивные и волновые потери в сбалансированном крейсерском полете при заданных коэффициенте подъемной силы, запасе продольной статической устойчивости и числе Маха. Форма в плане и распределение относительных толщин на крыле и оперении фиксированы. Методика состоит из прямого расчета компоновки крыло-оперение, основанного на решении уравнения для полного потенциала скорости, аналитического представления контуров сечений крыла и метода поиска экстремума. В качестве целевой функции выбрано произведение максимального локального числа Маха на поверхности крыла на индуктивное сопротивление компоновки крыло^оперение.
Ключевые слова: трансзвуковой крейсерский полет, компоновка крыло горизонтальное оперение, поиск экстремума.
1
1'2 2 1
1
2
Surface design of the wing with horizontal tail in transonic
cruise flight
The method of profiles design of the wing and horizontal tail configuration is considered. The minimum induced and wave losses in balanced transonic cruise flight for a fixed lift coefficient and Mach number are achieved. The wing and tail plan forms and thickness distribution are constant. The method consists of the direct aerodynamic calculation of the wing tail configuration based on the solution of a full velocity potential equation, analytical approximation of the wing sections profiles and the extremum seeking algorithm. The target function is the product of the maximum local Mach number and the induced drag of the wing tail configuration. The results of designing the wing tail configuration for a typical passenger aircraft are show.
Key words: transonic cruise flight, wing tail configuration, extremum seeking
Проблеме аэродинамического проектирования посвящено большое количество работ, простое перечисление которых займет не одну сотню страниц. Это говорит о том, что аэродинамическое проектирование является сложной и важной задачей, не имеющей одного определенного алгоритма решения. В монографии [1], в которой дан обзор методов аэродинамического проектирования, очень верно замечено, что «проектирование является конечной целью исследований, и все другие действия должны быть нацелены на ее достижение». Описанная в данной статье методика опирается на классическое представление проблемы проектирования [2, 3], в основе которой лежит прямой метод расчета аэродинамических
Ф) Теперин JI. Л., Теперина JI. Н., Мьё Тхейн, 2019
(с) Федеральное государственное автономное образовательное учреждение высшего образования
«Московский физико-технический институт (национальный исследовательский университет)», 2019
характеристик. Алгоритм аэродинамических) проектирования заключается в объединении нрямохх) метода расчета с двумя дополнительными методиками. Одна из них это параметризация формы поверхности, позволяющая изменять геометрическую форму летательших) аппарата, удовлетворяя принятым конструктивным ограничениям и гладкости вариаций. Вторая методика это а.;п'оритм поиска экстремума целевой функции, зависящей от неременных параметров задачи.
В качестве прямой задачи используется программа расчета трансзвуковохх) обтекания крыла, основанная на решении нелинейжнх) уравнения для полного потенциала скорости конечноразноетным методом [4|. Эта программа была модифицирована [5] таким образом, что появилась возможность рассчитывать не только изолированное крыло, но и крыло с расположенным вниз но потоку горизонтальным оперением. Приемлемая точность и высокая оперативность программы позволяют использовать ее в а.;п'оритме проектирования, при котором необходимо выполнять большое количество обращений к прямому расчету.
При заданной форме в плане аэродинамическая компоновка крыла определяется контурами сечений (профилями) и круткой. В отличие от обычного способа описания контуров сечений крыла табличной функцией в работе применено аналитическое представление вида [6|:
п
у(х) = -\f2px +
г=1
где р - радиус затупления передней кромки.
В работе [7] эта функция с полиномиальной частью 4-го порядка была использована для проектирования аэродинамической поверхности крыла ромбовидной формы в плане. В данной работе полиномиальная часть функции для описания отдельно верхних) и нижних) контуров сечений крыла увеличена до 7. На рис. 1 показана верхняя часть контура профиля с указанием трех точек, в которых выполняются условия для нахождения коэффициентов полинома а,г .
Рис. 1. Верхний контур профиля в сечении крыла
В точках контура 12 3 выполняются условия, указанные в таблице 1.
Для определения коэффициентов полинома в зависимости от заданных параметров контура необходимо решить систему обыкновенных линейных уравнений, в которой принято, что все параметры отнесены к хорде профиля с.
Таблица1
Условия для нахождения коэффициентов полинома
Номер точки Заданные параметры Обозначения Номер параметр
1 Радиус затупления Р 1
2 Максимум контура Ст 2
Положение максимума 3
Вторая производная Ут 4
3 Толщина задней кромки С1 5
Угол наклона задней кромки дг 6
Вторая производная У1 7
Площадь контура 5 8
аг ' 1 1 1 1 1 1 1 -1 С1 - /Др
(12 1 2 3 4 5 6 7 tan вг - л/р/2
аз 1 1хт 3т2 4т3 5Хт 6Хт 7хт р/2хт
а4 = /у2 ™3 хт Т6 7 хт С-т - л/2рХт
а5 0 2 6 12 20 30 42 у" + ^р/2/2
а6 0 2 бхт 12Хт 20хт 42Хт Ут + V р/(2хт)/2хт
а7 .1/2 1/3 1/4 1/5 1/6 1/7 1/8. Б - 2/3 * /2р
Теперь контур профиля можно описать не коэффициентами полинома, а естественными параметрами. Нижний контур представляется аналогичным образом. Поскольку сумма максимальных толщин верхнего и нижнего контуров задана (толщина «на просвет»), а толщину на задней кромке следует уравнять, чтобы не дублировать крутку сечения, то для описания всего профиля необходимо задать 15 параметров, включая крутку. Для описаний всей поверхности потребуется количество параметров, кратное числу сечений крыла. Для крыла, представленного тремя сечениями, количество параметров, подлежащих определению, становится равным 45. Аэродинамическая компоновка оперения принимается заданной и не варьируется, за исключением угла установки. Кроме того, в задаче появляются ограничения в виде заданного коэффициента подъемной силы и нулевого продольного момента компоновки (условие балансировки). От этих ограничений можно избавиться, модифицировав прямую задачу так, чтобы она решалась не при заданном угле атаки а и угле установки оперения а при заданном коэффициенте подъемной силы су и нулевом продольном моменте. Используя локальную линеаризацию аэродинамических характеристик относительно текущего угла атаки ао и текущего угла установки оперения §о, можно получить выражения для этих параметров, удовлетворяющих ограничениям по подъемной силе и продольному моменту:
5 = (су - Су0 + Сутхо/т")/(сёу - СуШ^/т"); а = -(тго + тгг)/т",
где су0 0 - коэффициент подъемной силы и продольного момента при а = ао и 5 = 5о соответственно, су,ту,су- производные коэффициента подъемной силы и продольного момента по углу атаки и углу отклонения оперения соответственно. Относительно большое количество искомых переменных, отсутствие ограничений и малое компьютерное время, потребное для вызова прямой задачи, позволяет применить наиболее простой способ поиска экстремума - покоординатный спуск. В этом методе каждая искомая переменная приобретает приращение, которое затем оценивается по величине целевой функции. На рис. 2 показаны четыре возможных варианта отклика целевой функции Р(хг) на приращение параметра х
На начальном этапе выбирается шаг приращения параметра равный 10% его значения по модулю. После обращения к прямой задаче с уменьшенным и увеличенным на шаг
значением параметра анализируется отклик целевой функции. Если отклик соответствует варианту 1 на рис. 2, то значение параметра увеличивается на выбранный шаг, если варианту 2 значение параметра уменьшается на выбранный шаг. Варианты 3 («локальный минимум») и 4 («седло») игнорируются, и данный параметр не изменяет своего значения. Вычисление отклика целевой функции повторяется для всех параметров до тех нор, пока не исчезнут варианты 1 и 2. После этого шаг уменьшается в два раза, и процесс повторяется с уменьшенным шагом. Процесс прекращается после того как приращение целевой функции не станет меньше заданного значения. Для исключения возможного попадания в локальный минимум следует использовать несколько начальных приближений для искомых параметров.
Рис. 2. Варианты отклика целевой функции на приращение значения параметра
Остается вопрос о выборе целевой функции. Наиболее простым будет выбор в качестве целевой функции сопротивления компоновки. Поскольку в данном исследовании силами вязкости ирснсбрсгастся, в качестве сопротивления выступают две компоненты: волновое сопротивление и индуктивное сопротивление. Сопротивление без разделения на эти компоненты можно определить интегрированием проекции сил давления на направление набегающей скорости. Однако практика расчетов показывает, что такой способ определения сопротивления не всех да обладает необходимой точностью и требует большого количества расчетных точек на поверхности компоновки. С другой стороны, индуктивное сопротивление хорошо определяется в плоскости Трефтца по распределению циркуляции. Если крыло и оперение находятся в одной плоскости, распределения циркуляции по этим элементам в плоскости Трефтца объединяются, и индуктивное сопротивление определяется с достаточной точностью по методу работы [8]. Что касается волнового сопротивления, то оно будет отсутствовать, если на поверхности компоновки не будет скачков уплотнения. Очевидно, чем меньше число Маха на поверхности крыла, тем менее вероятно возникновение скачков уплотнения. Поэтому в качестве целевой функции для снижения волнового сопротивления целесообразно использовать максимальное локальное число Маха.
Применим целевую функцию последовательно в виде максимального локального числа Маха и в виде волнового сопротивления к проектированию методом покоординатного спуска контура крыла очень большого удлинения (профиля) при заданной толщине 12%, числе Маха (0.77) и коэффициенте подъемной силы (0.5). На рис. 3 показан результат проектирования профиля с использованием в качестве целевой функции сопротивления и максимального локального числа Маха.
Рис. За. Проектирование контура профиля но двум целевым функциям
Рис. 36. Проектирование контура профиля по двум целевым функциям
Видно, что распределение коэффициента давления (рис. За) и форма контуров профилей (рис. 36) для обеих целевых функций отличаются незначительно. Однако коэффициент сопротивления меньше, когда он является целевой функцией и, соответственно, максимальное локальное число Маха, которое пропорционально коэффициенту давления, также меньше, если его выбирать в качестве целевой функции. У компоновки крыла конечного размаха
и оперения появляется еще одна компонента сопротивления индуктивное сопротивление. Чтобы объединить требование минимального волнового и индуктивного сопротивления, предлагается ввести целевую функцию в виде произведения максимального локального числа Маха на индуктивное сопротивление. Обе величины строго положительны и снижение любой из них будет приводить к снижению целевой функции. На наш взгляд, такой подход более эффективен, чем введение суммы с весовыми коэффициентами. Волновое сопротивление можно контролировать по распределению давления на поверхности крыла. Применим разработанную методику к модельному примеру компоновки крыла и горизонтального оперения, изображенной на рис. 4. Будем искать форму профилей в трех сечениях крыла: в плоскости симметрии, в изломе и на конце, а крутку в двух сечениях: в изломе и на конце. Будем определять также угол атаки и угол установки оперения вместе с остальными параметрами, обеспечивающими минимум локального числа Маха и индуктивного сопротивления при заданном коэффициенте подъемной силы (су = 0.5), числе Маха крейсерского сбалансированного полета (М 0.8) и диапазоне запаса продольной статической устойчивости (Axf = 00! )• Крутка крыла в плоскости симметрии фиксирована и равна 3°. Относительная толщина профилей крыла в опорных сечениях задана с = 0.14-0.10-0.09. У оперения, расположенного с крылом в одной плоскости и построенного по симметричным профилям МАСА-ООЮ, варьируется только угол установки ¿г.0..
) 30.2'
= 185.41 мг
САХ = 4.61 м
40,88 и
24,115 м
Шг*
5Г()= 43.55 м2
14.888 м
Рис. 4. Форма в плане и основные размеры компоновки крыло оперение
Для проектирования одного варианта аэродинамической компоновки крыла требуется не более 30 000 обращений к прямому расчету. Если учесть, что время одного расчета на современном персональном компьютере составляет 7 секунд, то процесс оптимизации занимает не более 58.33 часа. В табл. 2 приведены основные параметры, полученные в процессе оптимизации для различных запасов продольной статической устойчивости Axf.
Из табл. 2 следует, что поиск минимума целевой функции позволяет получить для любого запаса устойчивости небольшое приращение индуктивного сопротивления относительно крыла с эллиптическим распределением циркуляции:
Cximin = С2/ж\(\ = 9, 01).
Т а б .л и ц а 2
Параметры оптимизации
A% 0,05 0,15 0,25 0,35 0,45
M 0,8 0,8 0,8 0,8 0,8
Су 0,5 0,5 0,5 0,5 0,5
отрад 0,64 0,71 0,76 0,81 0,83
.оград 1 )53 1,30 0,99 0,84 0,51
Крутка в пл.сим. град. 3 3 3 3 3
Крутка в изломе град. -0,61 -0,54 -0,67 -0,67 -0,68
Крутка на конце град. -1,31 -1,26 -1,46 -1,54 -1,59
Су кр 0,5007 0,5108 0,5224 0,5329 0,5420
СуТ .0 -0,0007 -0,0108 -0,0224 -0,0329 -0,0420
™>z (су = 0) без г.о -0,1419 -0,1382 -0,1405 -0,1372 -0,1378
Cx^Cximin 1,006 1,005 1,005 1,006 1,006
Рис. 5. Распределение циркуляции по крылу и компоновке крыло оперение
На рис. 5 показан механизм снижения индуктивного сопротивления в зависимости от запаса устойчивости. При большом запасе устойчивости ( Axf = 0.45) для балансировки требуется значительная отрицательная подъемная сила оперения и соответствующее отрицательное распределение циркуляции, показанное на рис. 5 точечной линией. В процессе оптимизации на крыле образуется колоколообразное распределение циркуляции (толстая пунктирная линия), и суммарное распределение циркуляции крыла и оперения (черная сплошная линия) будет приближаться к эллиптическому распределению (тонкая пунктирная линия). Для балансировки при малом запасе устойчивости (Axf = 0.15) на оперении достаточно образовать совсем небольшую отрицательную подъемную силу. Соответствующее распределение циркуляции оперения показано на рис. 5 звездочками. Оптимальное решение для циркуляции крыла (черные кружки) учитывает, что сложение его циркуляции с циркуляцией оперения опять дает результат, близкий к эллипсу (полые
кружки). Представляет интерес сделать оценку индуктивного сопротивления компоновки, спроектированной для одного запаса устойчивости, при изменении этого запаса на другую величину.
Рис. 6. Приращение индуктивного сопротивления при переходе на другой запас
устойчивости
Нижняя сплошная кривая на рис. 6 означает приращение индуктивного сопротивления оптимизированной компоновки относительно минимально возможного, при заданном запасе устойчивости. Остальные кривые получены расчетом потерь индуктивших) сопротивления при замене запаса устойчивости относительно проектируемого до больших) или меньших) значения. При этом у компоновки изменяются только угол атаки и угол установки оперения, а форма профилей и крутка крыла не меняются. Например, пунктирная кривая, помеченная на графике как Axf = 0.45, получена расчетом индуктивного сопротивления компоновки, спроектированной для Axf = 0.45 при изменении центра тяжести таким образом, что запас устойчивости изменяется от 0.45 до 0.05. Проведенные расчеты показывают, что минимальные индуктивные потери при смене центровки возможны для вариантов со средним запасом устойчивости 0.25.
На рис. 7 показано распределение давления и контур среднего сечения крыла (с = 0.55), спроектированного для запаса устойчивости Axf = 0.45 (сплошные кривые) и Axf = 0.15 (пунктирные кривые). Результаты проектирования показывают, что при изменении запаса продольной статической устойчивости изменяются не только угол атаки, угол установки горизонтального оперения и углы крутки крыла (см. таблицу 2), но и форма профилей в сечениях крыла. То есть в перестройке распределения циркуляции вдоль размаха крыла участвуют все искомые параметры. На рис. 7 видно, что в средней части размаха профиль крыла имеет весьма тонкий хвостик. Это значит, что результат оптимизации может не соответствовать ограничениям, связанным с прочностью конструкции крыла. Чтобы утолстить хвостик, введем в опорном сечении в изломе крыла ограничение на угол наклона нижнего контура на задней кромке. Введем условие равенства этого параметра на верхнем и нижнем контурах. Проведя оптимизацию целевой функции при дополнительном ограничении, получим новое решение с утолщенным хвостиком.
На рис. 8 показано распределение коэффициента давления и контуров профилей (сплошная кривая) без ограничения на угол наклона нижнего контура на задней кромке и с ограничением (пунктирная кривая). Очевидно, что ограничение, снизив «подрезку»
Рис. 7. Раиределение коэффициента давления и контуров профилей в среднем сечении
крыла
Рис. 8. Распределение коэффициента давления и контура исходного профиля и утолщеного на задней кромке профиля
на нижнем контуре в изломе перераспределило нагрузку в корневые сечения крыла. В результате за утолщение хвостика придется «заплатить» увеличением индуктивного сопротивления на 0.6%.
Значительная стреловидность крыла позволяет рассмотреть увеличение числа Маха крейсерского полета до М 0.82. При увеличении числа Маха набегающего потока локальное число Маха на поверхности возрастает интенсивнее в концевых сечениях крыла.
Рис. 9. Распределение циркуляции вдоль размаха и коэффициента давления в корневом сечении крыла спроектированного для М 0,8 и М 0,82
По этой причине, чтобы уменьшить локальное число Маха, необходимо перераспределить нагрузку ближе к корню крыла. Это, в свою очередь, вызывает рост индуктивного сопротивления. В результате процесс оптимизации сходится к компромиссному решению. На рис. 9 показано распределение циркуляции вдоль размаха крыла и распределение коэффициента давления в сечении г = 0.1 для двух значений числа Маха набегающего потока. Видно, что повышение крейсерского числа Маха вызывает рост индуктивного сопротивления на 0.65% и появление скачков уплотнения небольшой интенсивности. Для более точной оценки влияния увеличенного числа Маха крейсерского полета на аэродинамические характеристики компоновки необходимо иметь надежную методику оценки волнового сопротивления.
Проведенное численное проектирование аэродинамической компоновки крыла и горизонтального оперения, находящегося с ним в одной плоскости, показало, что включение оперения в процесс оптимизации позволяет снизить индуктивные потери до 0.5%; относительно крыла с эллиптическим распределением циркуляции в широком диапазоне запаса продольной устойчивости. Это объясняется тем, что кол околообразное распределение циркуляции на крыле, которое позволяет ему преодолеть волновой кризис, компенсируется отрицательным распределением циркуляции на оперении. Сложение циркуляций крыла и оперения в процессе оптимизации приводит к снижению индуктивного сопротивления компоновки в целом. При изменении продольного запаса устойчивости проектирование производится заново, при этом изменяется угол атаки, угол установки горизонтального оперения, крутка крыла и форма профилей, образующих аэродинамическую поверхность. Расчеты показали, что если компоновку спроектировать с запасом устойчивости Axf = 0.25, то при сдвиге центра тяжести назад на 20%; С АХ, приращение индуктивного сопротивления будет не более 0.6%;.
Литература
1. Брутян М.А. Основы трансзвуковой аэродинамики. Москва : Наука, 2017. С. 175.
2. Болсуповский А.Л., Бузоверя Н.П., Барипов В.А., Скоморохов С.И., Чернавских Ю.Н. Проектирование аэродинамической компоновки пассажирских и транспортных самолетов // Полет, спец. вып. 2008. С. 16-23.
3. Болсуповский А.Л., Бузоверя Н.П., Кара,съ О.В., Ковалев В.Е. Развитие методов аэродинамического проектирования крейсеркой компоновки дозвуковых самолетов // Труды ЦАГИ. Вып. 2655.2002. С. 133-145.
4. Кара,съ О.В., Ковалёв В.Е. Применение обратного метода расчёта трёхмерного пограничного слоя к задаче обтекания крыла с учётом влияния вязкости // Учёные записки ЦАГИ. 1989. T. XX, № 5. С. 1 11.
5. Теперин Л.Л., Теперина Л.Н., Мьё Тхейн Методика расчета аэродинамических характеристик компоновки крыло-оперение в трансзвуковом крейсерском полете // Учёные записки ЦАГИ. 2018. T. XLIX, № 6. С. 28-34.
6. Кутищев Г.П., Теперин Л.Л. Применение аналитического представления контура профиля для решения обратной задачи аэродинамики // Численные методы аэродинамического проектирования // Труды ЦАГИ. 2003. Вып. 2655. С. 104-108.
7. Лазарев В.В., Павленко А.А., Разов А.А., Теперин Л.Л., Теперина Л.Н. Аэродинамическое проектирование летательного аппарата ромбовидной формы в плане // Ученые записки ЦАГИ. 2011. T. X I'll. № 4. С. 30-37.
8. Теперин Л.Л., Притуло Т.М., Орфинежад Ф.Э., Мьё Тхейн Средства снижения индуктивного сопротивления // Труды МФТИ. 2017. Т. 9, № 4. С. 94-105.
References
1. Brutyan M.A. Basics of transonic aerodynamics. Moscow : Science. 2017. P. 175.
2. Bolsunovsky A.L., Buzoverya N.P., Barinov B.A., Skomorokh S.I., Chernavsky U.N. Design of aerodynamic configuration of passenger and transport aircraft. Flight, spec, issue. 2008. P. 16-23.
3. Bolsunovsky A.L., Buzoverya N.P., Karas O.B., Kovalev B.E. Development of methods of aerodynamic design of subsonic aircraft cruising layout. Trudy TsAGi. I. 2655.2002. P. 133145.
4. Kara O.B., Kovalev B.E. Application of the inverse method of calculation of the three-dimensional boundary layer to the problem of flow around the wing taking into account the influence of viscosity. Science Notes TsAGI. 1989. T. XX, N 5. P. 1-11.
5. Teperin L.L., Teperina L.N., Myo Thyen Method of calculation of aerodynamic characteristics of wing-tail layout in transonic cruising flight. Science Notes TsAGI. 2018. T. XLIX, N 6. P. 28-34.
6. Kutiev O.P., Teperin L.L. Application of analytical representation of the profile contour for solving the inverse problem of aerodynamics. In Numerical methods of aerodynamic design. Trudy TsAGI. 2003. I. 2655. P. 104-108.
7. Lazarev B.B., Pavlenko A.A., Razov A.A., Teperin L.L., Teperina L.N. Aerodynamic design of a diamond-shaped aircraft in terms of. Science Notes TsAGI. 2011. T. XTII, N 4. P. 30-37.
8. Teperin L.L., Pritulo T.M., Orfininezhad F.A., Myo Thyen Resources of inductive drag diminution. Proceedings of MIPT. 2017. V 9, N 4. P. 94-105.
Поступим в редакцию 14-02.2019