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

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

• (6) Теперин Л. Л., Теперина Л. Н., Мьё Тхейн, 2019

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

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

При заданной форме в плане аэродинамическая компоновка крыла определяется контурами сечений (профилями) и круткой. В отличие от обычного способа описания контуров сечений крыла табличной функцией в работе применено аналитическое представление вида. [6]:

п у(х) = у/ 2рх + У^ а_гх^г, г=1

где р - радиус затупления передней кромки.

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

Рис. 1. Верхний контур профиля в сечении крыла.

В точках контура. 1-2-3 выполняются условия, указанные в таблице 1.

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

Таблица!

Условия для нахождения коэффициентов полинома

Номер точки	Заданные параметры	Обозначения	Номер параметр
1	Радиус затупления	Р	1
2	Максимум контура	с т	2
	Положение максимума	ж т	3
	Вторая производная	” У т	4
3	Толщина задней кромки	сі	5
	Угол наклона задней кромки	Өі	6
	Вторая производная	У1	7
	Площадь контура	S	8

01		⎡1	1	1	1	1	1	1⎤	- 1	сі - V2P
02		1	2	3	4	5	6	7		tan Ө1 — р pf2
03		1	2ж т	2 ° жт	3 т	4 и ж т	бжт	6 1 ж т		— V Р/2жт
04	⎥=	ж т	2 ж т	ж т	4 ж т	5 ж т	жт	7 жт		с т  V2pжт
05		0	2	6	12	20	30	42		у1 + Vp/2/2
06		0	2	6ж т	12жт	20жт	зожт	42жт		Ут + //р/(2жт)/2жт
О7		[ 1/2	1/3	1/4	1/5	1 / 6	1/7	1/8 J		S — 2/3 * V2P

Теперь контур профиля можно описать не коэффициентами полинома, а естественными параметрами. Нижний контур представляется аналогичным образом. Поскольку сумма максимальных толщин верхнего и нижнего контуров задана (толщина «на просвет»), а толщину на задней кромке следует уравнять, чтобы не дублировать крутку сечения, то для описания всего профиля необходимо задать 15 параметров, включая крутку. Для описаний всей поверхности потребуется количество параметров, кратное числу сечений крыла. Для крыла, представленного тремя сечениями, количество параметров, подлежащих определению, становится равным 45. Аэродинамическая компоновка оперения принимается заданной и не варьируется, за исключением угла установки. Кроме того, в задаче появляются ограничения в виде заданного коэффициента подъемной силы и нулевого продольного момента компоновки (условие балансировки). От этих ограничений можно избавиться, модифицировав прямую задачу так, чтобы она решалась не при заданном угле атаки а и угле установки оперения 6, а при заданном коэффициенте подъемной силы с у и нулевом продольном моменте. Используя локальную линеаризацию аэродинамических характеристик относительно текущего угла атаки ао и текущего угла установки оперения 6о, можно получить выражения для этих параметров, удовлетворяющих ограничениям по подъемной силе и продольному моменту:

6 = (с_у - с_у 0 + с “ т_го/т “ )/(с⁵_у — с “ т ® /т^); а = —(т_го + т ® )/т^ где с у 0, т_г0 - коэффициент подъемной силы и продольного момента при а = ао и 6 = 6о соответственно, с“,т“,сУ ,т ® - производные коэффициента подъемной силы и продольного момента по углу атаки и углу отклонения оперения соответственно. Относительно большое количество искомых переменных, отсутствие ограничений и малое компьютерное время, потребное для вызова прямой задачи, позволяет применить наиболее простой способ поиска экстремума - покоординатный спуск. В этом методе каждая искомая переменная приобретает приращение, которое затем оценивается по величине целевой функции. На рис. 2 показаны четыре возможных варианта отклика целевой функции F (ж^) на приращение параметра жъ

На начальном этапе выбирается шаг приращения параметра равный 10% его значения по модулю. После обращения к прямой задаче с уменьшенным и увеличенным на шаг значением параметра анализируется отклик целевой функции. Если отклик соответствует варианту 1 на рис. 2, то значение параметра увеличивается на выбранный шаг, если варианту 2 - значение параметра уменьшается на выбранный шаг. Варианты 3 («локальный минимум») и 4 («седло») игнорируются, и данный параметр не изменяет своего значения. Вычисление отклика целевой функции повторяется для всех параметров до тех пор, пока не исчезнут варианты 1 и 2. После этого шаг уменьшается в два раза, и процесс повторяется с уменьшенным шагом. Процесс прекращается после того как приращение целевой функции не станет меньше заданного значения. Для исключения возможного попадания в локальный минимум следует использовать несколько начальных приближений для искомых параметров.

Рис. 2. Варианты отклика целевой функции на приращение значения параметра

Остается вопрос о выборе целевой функции. Наиболее простым будет выбор в качестве целевой функции сопротивления компоновки. Поскольку в данном исследовании силами вязкости пренебрегается, в качестве сопротивления выступают две компоненты: волновое сопротивление и индуктивное сопротивление. Сопротивление без разделения на эти компоненты можно определить интегрированием проекции сил давления на направление набегающей скорости. Однако практика расчетов показывает, что такой способ определения сопротивления не всегда обладает необходимой точностью и требует большого количества расчетных точек на поверхности компоновки. С другой стороны, индуктивное сопротивление хорошо определяется в плоскости Трефтца по распределению циркуляции. Если крыло и оперение находятся в одной плоскости, распределения циркуляции по этим элементам в плоскости Трефтца объединяются, и индуктивное сопротивление определяется с достаточной точностью по методу работы [8]. Что касается волнового сопротивления, то оно будет отсутствовать, если на поверхности компоновки не будет скачков уплотнения. Очевидно, чем меньше число Маха на поверхности крыла, тем менее вероятно возникновение скачков уплотнения. Поэтому в качестве целевой функции для снижения волнового сопротивления целесообразно использовать максимальное локальное число Маха.

Применим целевую функцию последовательно в виде максимального локального числа Маха и в виде волнового сопротивления к проектированию методом покоординатного спуска контура крыла очень большого удлинения (профиля) при заданной толщине 12%, числе Маха (0.77) и коэффициенте подъемной силы (0.5). На рис. 3 показан результат проектирования профиля с использованием в качестве целевой функции сопротивления и максимального локального числа Маха.

Рис. За. Проектирование контура профиля по двум целевым функциям

Рис. 36. Проектирование контура профиля по двум целевым функциям

Видно, что распределение коэффициента давления (рис. За) и форма контуров профилей (рис. 36) для обеих целевых функций отличаются незначительно. Однако коэффициент сопротивления меньше, когда он является целевой функцией и, соответственно, максимальное локальное число Маха, которое пропорционально коэффициенту давления, также меньше, если его выбирать в качестве целевой функции. У компоновки крыла конечного размаха и оперения появляется еще одна компонента сопротивления - индуктивное сопротивление. Чтобы объединить требование минимального волнового и индуктивного сопротивления, предлагается ввести целевую функцию в виде произведения максимального локального числа Маха на индуктивное сопротивление. Обе величины строго положительны и снижение любой из них будет приводить к снижению целевой функции. На наш взгляд, такой подход более эффективен, чем введение суммы с весовыми коэффициентами. Волновое сопротивление можно контролировать по распределению давления на поверхности крыла. Применим разработанную методику к модельному примеру компоновки крыла и горизонтального оперения, изображенной на рис. 4. Будем искать форму профилей в трех сечениях крыла: в плоскости симметрии, в изломе и на конце, а крутку - в двух сечениях: в изломе и на конце. Будем определять также угол атаки и угол установки оперения вместе с остальными параметрами, обеспечивающими минимум локального числа Маха и индуктивного сопротивления при заданном коэффициенте подъемной силы (су = 0.5), числе Маха крейсерского сбалансированного полета (М = 0.8) и диапазоне запаса продольной статической устойчивости (Ажу = 0^45)• Крутка крыла в плоскости симметрии фиксирована и равна 3°. Относительная толщина профилей крыла в опорных сечениях задана с = 0.14-0.10-0.09. У оперения, расположенного с крылом в одной плоскости и построенного по симметричным профилям NACA-0010, варьируется только угол установки ф.0..

= 185.41 м

САХ = 4.61 м

737.2°

Sro= 43.55 м2

14.888 м

Рис. 4. Форма, в плане и основные размеры компоновки крыло-оперение

Для проектирования одного варианта, аэродинамической компоновки крыла, требуется не более 30 000 обращений к прямому расчету. Если учесть, что время одного расчета, на. современном персональном компьютере составляет 7 секунд, то процесс оптимизации занимает не более 58.33 часа. В табл. 2 приведены основные параметры, полученные в процессе оптимизации для различных запасов продольной статической устойчивости Ажу.

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

C ximin — С у /хХ^Х — 9, 01).

Т а б л и ц а 2

Параметры оптимизации

А	0,05	0,15	0,25	0,35	0,45
м	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
5Г .оград	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
mz(су — 0) без г.о	-0,1419	-0,1382	-0,1405	-0,1372	-0,1378
C x ^(C ximin	1,006	1,005	1,005	1,006	1,006

Рис. 5. Распределение циркуляции по крылу и компоновке крыло-оперение

На рис. 5 показан механизм снижения индуктивного сопротивления в зависимости от запаса устойчивости. При большом запасе устойчивости ( Ажу — 0.45) для балансировки требуется значительная отрицательная подъемная сила оперения и соответствующее отрицательное распределение циркуляции, показанное на рис. 5 точечной линией. В процессе оптимизации на крыле образуется колоколообразное распределение циркуляции (толстая пунктирная линия), и суммарное распределение циркуляции крыла и оперения (черная сплошная линия) будет приближаться к эллиптическому распределению (тонкая пунктирная линия). Для балансировки при малом запасе устойчивости (Ажу — 0.15) на оперении достаточно образовать совсем небольшую отрицательную подъемную силу. Соответствующее распределение циркуляции оперения показано на рис. 5 звездочками. Оптимальное решение для циркуляции крыла (черные кружки) учитывает, что сложение его циркуляции с циркуляцией оперения опять дает результат, близкий к эллипсу (полые кружки). Представляет интерес сделать оценку индуктивного сопротивления компоновки, спроектированной для одного запаса устойчивости, при изменении этого запаса на другую величину.

Рис. 6. Приращение индуктивного сопротивления при переходе на другой запас устойчивости

Нижняя сплошная кривая на рис. 6 означает приращение индуктивного сопротивления оптимизированной компоновки относительно минимально возможного, при заданном запасе устойчивости. Остальные кривые получены расчетом потерь индуктивного сопротивления при замене запаса устойчивости относительно проектируемого до большего или меньшего значения. При этом у компоновки изменяются только угол атаки и угол установки оперения, а форма профилей и крутка крыла не меняются. Например, пунктирная кривая, помеченная на графике как Ажу = 0.45, получена расчетом индуктивного сопротивления компоновки, спроектированной для Ажу = 0.45 при изменении центра тяжести таким образом, что запас устойчивости изменяется от 0.45 до 0.05. Проведенные расчеты показывают, что минимальные индуктивные потери при смене центровки возможны для вариантов со средним запасом устойчивости 0.25.

На рис. 7 показано распределение давления и контур среднего сечения крыла (г = 0.55), спроектированного для запаса устойчивости Ажу = 0.45 (сплошные кривые) и Ажу = 0.15 (пунктирные кривые). Результаты проектирования показывают, что при изменении запаса продольной статической устойчивости изменяются не только угол атаки, угол установки горизонтального оперения и углы крутки крыла (см. таблицу 2), но и форма профилей в сечениях крыла. То есть в перестройке распределения циркуляции вдоль размаха крыла участвуют все искомые параметры. На рис. 7 видно, что в средней части размаха профиль крыла имеет весьма тонкий хвостик. Это значит, что результат оптимизации может не соответствовать ограничениям, связанным с прочностью конструкции крыла. Чтобы утолстить хвостик, введем в опорном сечении в изломе крыла ограничение на угол наклона нижнего контура на задней кромке. Введем условие равенства этого параметра на верхнем и нижнем контурах. Проведя оптимизацию целевой функции при дополнительном ограничении, получим новое решение с утолщенным хвостиком.

На рис. 8 показано распределение коэффициента давления и контуров профилей (сплошная кривая) без ограничения на угол наклона нижнего контура на задней кромке и с ограничением (пунктирная кривая). Очевидно, что ограничение, снизив «подрезку»

Рис. 7. Рапределение коэффициента давления и контуров профилей в среднем сечении крыла

Рис. 8. Распределение коэффициента давления и контура исходного профиля и утолщеного на задней кромке профиля на нижнем контуре в изломе перераспределило нагрузку в корневые сечения крыла. В результате за утолщение хвостика придется «заплатить» увеличением индуктивного сопротивления на 0.6%.

Значительная стреловидность крыла позволяет рассмотреть увеличение числа Маха крейсерского полета до М = 0.82. При увеличении числа Маха набегающего потока локальное число Маха на поверхности возрастает интенсивнее в концевых сечениях крыла.

Рис. 9. Распределение циркуляции вдоль размаха и коэффициента давления в корневом сечении крыла спроектированного для М = 0,8 и М = 0,82

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

Проведенное численное проектирование аэродинамической компоновки крыла и горизонтального оперения, находящегося с ним в одной плоскости, показало, что включение оперения в процесс оптимизации позволяет снизить индуктивные потери до 0.5% относительно крыла с эллиптическим распределением циркуляции в широком диапазоне запаса продольной устойчивости. Это объясняется тем, что колоколообразное распределение циркуляции на крыле, которое позволяет ему преодолеть волновой кризис, компенсируется отрицательным распределением циркуляции на оперении. Сложение циркуляций крыла и оперения в процессе оптимизации приводит к снижению индуктивного сопротивления компоновки в целом. При изменении продольного запаса устойчивости проектирование производится заново, при этом изменяется угол атаки, угол установки горизонтального оперения, крутка крыла и форма профилей, образующих аэродинамическую поверхность. Расчеты показали, что если компоновку спроектировать с запасом устойчивости Ажу = 0.25, то при сдвиге центра тяжести назад на 20% САХ, приращение индуктивного сопротивления будет не более 0.6%.