Определение аэродинамических характеристик беспилотного летательного аппарата самолетного типа аналитическими методами
Приводится методика расчета аэродинамических характеристик беспилотного летательного аппарата (БПЛА) самолетного типа с помощью аналитических методов. Методика позволяет определять зависимости аэродинамических характеристик на рабочих углах атаки по известным моделям, а также экстраполировать результаты численных расчетов или экспериментов на весь диапазон рабочих углов атаки. Использование данной методики позволит форсировать этапы концептуального и предварительного проектирования БПЛА, определить режимы полета с максимальными значениями аэродинамических критериев эффективности. С использованием данной методики спроектированы БПЛА и проведено сравнение результатов расчета с результатами трехмерного моделирования аэродинамики планеров БПЛА в программных пакетах OpenFOAM и ANSYS Fluent на ресурсах суперкомпьютера Томского государственного университета СКИФ Cyberia.
Determination of aerodynamic characterisitcs of fixed-wing unmanned aerial vehicle by analytical techniques.pdf В данной работе рассматривается методика расчета аэродинамических характеристик летательного аппарата самолетного типа, применимая к беспилотным летательным аппаратам (БПЛА). Аэродинамические характеристики планера БПЛА в совокупности с его взлетной массой и характеристиками его силовой установки задают оптимальные режимы полета, при которых достигаются наибольшие или наименьшие значения целевых функций: дальности, продолжительности полета, вертикальной скорости, набора высоты и других величин, определяющих экономическую целесообразность эксплуатации. В частности, общими критериями целевых функций являются аэродинамическое качество K и критерий , на основе которых в работе [1] показано определение оптимальных режимов крейсерского полета летательного аппарата. Поэтому задача определения аэродинамических характеристик БПЛА еще на этапе предэскизного и эскизного проектирования имеет важное значение для выбора наиболее эффективных конфигураций летательного аппарата. Для решения такой задачи широко применяются методы оценок, численные исследования аэродинамики взамен дорогостоящих экспериментальных методов исследования. Инженерные методы оценок аэродинамических характеристик основаны на классической теории крыла, общих моделях аэродинамики некоторых тел, автомодельных решениях и наиболее подробно описываются в работах [2-4]. Численные решения задачи основываются на известных подходах панельных методов, методов вычислительной гидромеханики и положениях аэродинамики летательных аппаратов. На основе методов вычислительной гидромеханики реализованы методы оптимального аэродинамического проекти-113 Механика / Mechanics рования, использующие генетические алгоритмы [5, 6]. Для всех этих методов требуется геометрическая модель планера, параметры которого на этапе предварительного проектирования неизвестны. Ввиду трудоемкости прямого решения, а иногда и невозможности определения геометрических параметров при требуемой аэродинамике, определение этих параметров планера производится методом последовательных приближений, при котором количество циклов эскизного проектирования от начального приближения до получения приближения летательного аппарата, удовлетворяющего эксплуатационным требованиям, может быть велико. Существующие и доступные методы расчета первого приближения БПЛА основаны на методах проектирования пилотируемых летательных аппаратов и не учитывают различий конструкции, энергетических и силовых установок, оборудования, характерных чисел Рейнольдса и других коэффициентов. Поэтому целью данной работы является разработка методов для определения аэродинамических характеристик, оптимальных режимов полета без привязки к геометрическим моделям планера БПЛА, а лишь к наиболее общим геометрическим характеристикам. Геометрическое моделирование должно быть следующим за предварительным расчетом этапом проектирования, и его параметры должны определяться на основе требуемых аэродинамических характеристик. Постановка и метод решения задачи Рассматривается задача полета БПЛА, включающего основные фазы: взлет, набор высоты, горизонтальный полет, снижение, посадку. При проведении расчетов предполагалось, что фазы полета независимы друг от друга и представляют собой процессы равномерного установившегося движения, за исключением фазы взлета. Поэтому переходными процессами в рамках одной фазы полета пренебре-гается. Между различными фазами полета переходные процессы учитываются при расчете энергии бортового источника. Основные силы, воздействующие на БПЛА в полете, приведены на рис. 1. Будем рассматривать рабочий диапазон углов атаки, при которых угол атаки а считается малым. Динамика бокового движения не рассматривается в силу малости аэродинамических боковых сил по сравнению с продольными и поперечными силами. Также исключается из рассмотрения динамика вращательного движения, так как моменты, воздействующие на планер БПЛА, в равновесном состоянии взаимно сбалансированы. Для удобства представим аэродинамические силы через продольную и нормальную силы связанной системы координат, X и Y соответственно: X = Ya а - Xa; (1) Y = Y а + Ха а. (2) Тогда уравнения баланса сил в проекциях на связанную систему координат имеют вид: P +X - G sinu = 0; (3) Y - G cos u = 0. (4) В фазе взлета силовая установка БПЛА совершает работу, направленную на увеличение скорости до Vmin и набор до минимальной безопасной высоты Hmin. С учетом КПД во взлетном режиме цс.у.взл. энергию для взлета приближенно можно записать: 114 Исмаилов К.К. Определение аэродинамических характеристик W„, V2 gH Ш1П + -ү- (5) Рис. 1. Системы координат и силы, действующие на БПЛА. Ox&Vg - земная система координат, Oxaya - скоростная система координат, Оху - связанная система координат, а - угол атаки, и - угол тангажа, V - воздушная скорость, Ха - сила аэродинамического сопротивления, Ya - аэродинамическая подъемная сила, G - вес БПЛА, Р - сила тяги силовой установки Fig. 1. Coordinate systems and the forces acting on a UAV. Ox&Vg is the Earth fixed coordinate system, Oxaya is the wind-axes coordinate system, Oxy is the aircraft body fixed coordinate system, a is the angle of attack, n is the pitch angle, V is the airspeed, Xa is the aerodynamic drag force, Ya is the aerodynamic lift force, G is the UAV weight, and P is the force of propulsion Длительность фазы набора /паб. определяется заданной высотой H, углом набора и - а и скоростью Ѵаб.: _ AH . наб. -наб. (sin и - acos и) ’ AH = H - Hmin - H , где Hi - высота точки взлета БПЛА над уровнем моря. Тогда работа силовой установки с учетом КПД в режиме набора ^с.у.наб. W с.у.наб. 1 1 c cos и 1 acosu GAH, П с.у.наб. С Sinu sinu (6) где Cx, Cy - коэффициенты сил X, Y. Работа силовой установки в фазе горизонтального полета, продолжительностью /кр. при скорости Ѵкр., с учетом КПД силовой установки рс.у.кр. 1 С - (7) W„, Пс x-GV t кр. кр. С учетом (i)-(2) и вводя аэродинамическое качество как K = суа/сш выражения (6)-(7) примут вид: 115 Механика / Mechanics W с.у.наб. - 1 1 1+■ Пс.у.наб. W с.у.кр. _ 1 f1 Пс.у.кр. 1К 1 cos и + asinu K sinu - acosu GAH кр. кр (8) (9) Из выражений (8)-(9) видно, что при фиксированных и, V и tKp. стигается при условии максимума аэродинамического качества K. через Y, равный G в горизонтальном полете, перепишем выражение 1 с„ „ ч G15 минимум доВыражая Ѵкр. (7): п с15 Іс.у.кр. ya (1 - Ka) V pS0 t„ где р - плотность воздуха, So - базовая площадь крыла. В последнем выражении достигается минимальное значение мощности горизонтального полета ^с.у.кр. при условии максимума критерия Kjd = d /cxa . Предполагается, что при снижении БПЛА переходит в режим планирования, поэтому работа силовой установки равна нулю. Однако при проведении расчетов необходимо выполнить учет потребления целевой нагрузки и бортового оборудования. Считая их мощность потребления постоянной, для расчета потребной энергии достаточно вычислить продолжительность снижения: AH =- V ( план. V sinu - acosu У AH = H - H , где Ѵішан. - скорость планирования, Н2 - высота точки приземления над уровнем моря. Из выражений (1)-(4) для фазы снижения с учетом коэффициентов сил получим 1 - a = tg(-u). Суммарное энергопотребление силовой установки БПЛА в полете (5), (8) и (9) минимально при максимуме параметра K, требуемая мощность минимальна при максимуме критерия dC, что при фиксированном запасе энергии соответствует наибольшей дальности в первом случае и максимальной продолжительности полета во втором. Эти критерии определяют совершенство аэродинамической компоновки планера БПЛА и являются показателями его эффективности. В работе [7] отмечается, что индуктивное сопротивление крыла конечного размаха минимально при эллиптическом распределении циркуляции по размаху. Для такого крыла справедливы выражения, определяющие коэффициент подъемной силы cya и коэффициент сопротивления cxa: (10) (ii) Сya = СУ (a - a0 ); 2 с с = с +-- . xa xp п nek Далее в тексте при обозначении аэродинамических сил и их коэффициентов индекс a будем опускать, т.е. символами X, Y, сх, Су будем обозначать аэродина-116 Исмаилов К.К. Определение аэродинамических характеристик мические силы и их коэффициенты вместо сил связанной системы координат. В (10) cay определялось в соответствии с рекомендациями работ [3, 4]: < = 2п fl + 2 У Угол атаки а0 - угол, при котором подъемная сила равна нулю (нулевая подъемная сила): а0Ізв = а0І2Б ; “o|2D =-1.04/ , где а0|зс - угол атаки нулевой подъемной силы крыла конечного размаха, a0|2D -угол атаки нулевой подъемной силы в плоской задаче обтекания профиля крыла при характерных для крыла БПЛА числах Рейнольдса (Re = 105 + 2 -106), f -относительная кривизна профиля в процентах. В выражении (13) схр - коэффициент профильного сопротивления планера БПЛА, X - удлинение крыла, e - коэффициент Освальда, отражающий интерференцию крыла с фюзеляжем, мотогондолами и другими элементами на крыле. Коэффициент профильного сопротивления в первом приближении можно представить в виде суммы сопротивлений отдельных элементов планера: с = Ѵ с . (12) xp хр.элем ѵ / Коэффициент профильного сопротивления для крыла Схркрыло и оперения Схропер берется из поляр аэродинамических профилей. Для фюзеляжа, мотогондол и других выступающих элементов коэффициент профильного сопротивления берется из табличных данных. В общем случае коэффициент состоит из сопротивления давления и сопротивления трения: с хр.фюз. S„ = с -- + с -- хр.фюз.я ^ xp.фюз.Т S0 S0 Коэффициент сопротивления давления c зависит от формы атакующих поверхностей с относительной площадью Sa, коэффициент сопротивления трения с^фюзт более консервативен и зависит от степени шероховатости контактируемой (омываемой) поверхности Sп и локальных чисел Рейнольдса Re/. Форма и размеры фюзеляжа и мотогондол задаются формой, расположением и габаритами целевой нагрузки, двигателей, бортового источника питания и другого оборудования. Тем не менее в хорошо спроектированном планере с развитым аэродинамически эффективным корпусом сопротивление давления невелико и в первом приближении может быть учтено только для выступающих элементов целевой нагрузки, антенн и т.д. Сопротивление трения может быть найдено на основе сопротивления плоской пластины с использованием номограмм (рис. 2) или известных зависимостей для гладких и шероховатых пластин [8]: с f 0.074 Re 0 2 (5 -105 < Re < 107) ; (13) Mb '~S~ 117 Механика / Mechanics Re, = H. V где характерная ширина фюзеляжа АІ и его длина b, ѵ - вязкость воздуха. Для корпусов БПЛА, изготовленных из композиционных материалов методом формования в матрицах, относительная шероховатость равна b/ks ~ ІО4 -гІО6. Рис. 2. Закон сопротивления плоской пластины с песочной шероховатостью; полный коэффициент трения (/ = b , U^= V, для поверхностей БПЛА l/ks - 104 н- 10б) [8. С. 587] Fig. 2. Resistance law for a sand-roughened plate; a coefficient of total skin friction (/ = b , Ur, = V , for UAVs’ surfaces l/ks ~ ІО4 н-106) [8. P. 653] Интерференция крыла с элементами, расположенными или занимающими часть площади крыла, приводит к неравномерности скоса потока за крылом, следовательно, к увеличению индуктивного сопротивления. Данный эффект учитывается в коэффициенте Освальда, который определяется следующим образом: e = S„ So + Е К: A,Jf (14) где дополнительные слагаемые рядом с базовой площадью крыла S отражают площадь частей крыла, занятых выступающими на поверхности крыла элементами, и площадь поверхности крыла ниже по потоку от этих элементов. Коэффици-118 Исмаилов К.К. Определение аэродинамических характеристик ент расположения элементов относительно длины местной хорды кы и участок размаха крыла Д/ в первом приближении считаются консервативными по отно шению к площади и удлинению крыла. Тогда аэродинамическое качество планера с учетом (11) будет иметь вид: * =-L c c + x nek (15) Фиксируя величину подъемной силы Y и дифференцируя выражение (15) по воздушной скорости V, получим, что при некотором значении коэффициента подъемной силы c * достигается максимум К * = Ктах, которому соответствует воздушная скорость V *: сУ = VcxPnek; (16) * 1 Inek К = - 2 \\І с (17) (18) , * „ (19) V су PS0 Выражения (16)-(19) соответствуют режиму полета с наибольшим аэродинамическим качеством К = Кшах . Аналогичным образом, выполняя оценку критерия эффективности К^С, получим: с/* =V 3cxpnek; (20) ( ka/s ) = (nek /2; V * = 27 j^-** "^3 I nek 4 \\ c- («{Г,)'' = ^ (nek )32; V = 27 1 „ = ~1= V и 0.76V . (21) (22) (23) jc/*pSo ^3 Выражения (20)-(23) соответствуют режиму полета с наибольшим значением критерия Kfc =( КҒ, )_• Верификация метода и результаты расчетов С использованием приведенной методики спроектированы модели БПЛА Прото-тип-2Э и Прототип-2Т. Обе модели БПЛА выполнены в виде аэродинамической компоновки «летающее крыло» и должны выполнять полет в течение 120 мин и нести полезную нагрузку массой 2 кг с потребляемой мощностью 60 Вт. 119 Механика / Mechanics Геометрические характеристики БПЛА определены в соответствии с [9] и приведены в табл. 1. Таблица 1 Характеристики БПЛА Прототип-2Э и Прототип-2Т БПЛА Пр-2Э Пр-2Т Масса взлетная, кг 7.3 7.2 Площадь крыла, м2 0.663 0.693 Размах крыла, м 2 2 САХ, м 0.33 0.35 Удлинение крыла 5.91 5.94 Коэффициент Освальда 0.86 0.76 Положение ЦТ, % САХ 15 15 Минимальная скорость теоретическая, м/с 14 14 Скорость наибольшего качества (оценка), м/с 22 20 Аэродинамическое качество K (оценка) 16.7 16.0 Критерий эффективности (оценка) 12.6 13.6 Результаты эскизного проектирования БПЛА Прототип-2Э представлены в работе [10]. Крыло летательного аппарата обладало неустойчивостью в срыве и требовало пересмотра конструкции. При дальнейшей работе над эскизной моделью БПЛА Прототип-2Э из соображений безопасности полетов были внесены изменения в конструкцию крыла: сменились контур крыла, расположение оперения. По результатам работы была спроектирована рабочая модель БПЛА Прототип-2Т. Для проверки спроектированных моделей БПЛА проводилось трехмерное исследование аэродинамики. Скорости полета БПЛА соответствуют числам Рейнольдса 2 -105 + 9 -105. Рассматриваемое течение описывается уравнениями Рейнольдса [11], для замыкания которых использовалась модель турбулентности SST &-ю [12]. Для получения единственного решения на границах расчетной области ставились граничные условия в соответствии с рекомендациями [12]. Решение задачи проводилось на вычислительных узлах суперкомпьютера СКИФ Cyberia в программных пакетах ANSYS Fluent и OpenFOAM [13]. Для уравнений переноса субстанции в обоих решателях использовались алгоритм SIMPLE, описанный в работе [14], и численные схемы второго порядка аппроксимации для производных. В решателе ANSYS Fluent использовалась тетраэдральная расчетная сетка, полученная с помощью модуля ANSYS Meshing, в OpenFOAM - гекса-эдральная расчетная сетка, полученная с помощью утилиты SnappyHexMesh. На рис. 3, 4 приведены зависимости коэффициентов подъемной силы и силы сопротивления, аэродинамического качества и критерия . В рабочем диапазоне углов атаки планера БПЛА наблюдается хорошее совпадение результатов расчетов. С ростом углов атаки наблюдается отклонение расчета коэффициента сопротивления от результатов моделирования. Вероятно, это связано с тем, что распределение подъемной силы на больших углах атаки отличается от эллиптического закона, и индуктивное сопротивление крыльев БПЛА возрастает быстрее, чем у крыла с эллиптическим распределением. Также наблюдается ошибка определения положения угла атаки наибольшего аэродинамического качества K в 1-2° и максимума критерия - в 3-4°. Ошибки результа тов расчетов приведены табл. 2, 3. 120 Исмаилов К.К. Определение аэродинамических характеристик Рис. 3. Поляра и точки максимальных значений целевых функций для БПЛА Прототип-2Э + + + - данная методика, ▲▲▲ - Fluent, о о о - OpenFOAM Fig. 3. Polar and maximum of efficiency criteria for UAV Prototype-2E + + + - the proposed approach, ▲▲▲ - Fluent, and о о о - OpenFOAM Рис. 4. Поляра и точки максимальных значений целевых функций для БПЛА Прототип-2Т + + + - данная методика, ▲▲▲ - Fluent, о о о - OpenFOAM Fig. 4. Polar and maximum of efficiency criteria for UAV Prototype-2T + + + - the proposed approach, ▲▲▲ - Fluent, and о о о - OpenFOAM 121 Механика / Mechanics Т аблица 2 Характеристики БПЛА Прототип-2Э Параметр Су K КҒ> Среднеквадратичное отклонение 0.001 0.003 0.369 0.242 Ошибка при Vmin, % 5.0 0.3 5.2 5.4 Ошибка при V, % 1.9 1.6 0.3 0.5 Ошибка при V*, % 1.5 1.1 0.2 0.3 Таблица 3 Характеристики БПЛА Прототип-2Т Параметр С* Су K КҒ> Среднеквадратичное отклонение 0.002 0.009 0.526 0.588 Ошибка при Vmin, % 8.9 3.9 14.1 16.3 Ошибка при V*, % 0.2 2.1 2.2 3.3 Ошибка при V*, % 4.9 2.8 8.1 9.6 Заключение Приведено теоретическое исследование метода расчета аэродинамических характеристик беспилотного летательного аппарата. Методика позволяет установить основные аэродинамические характеристики БПЛА и при решении уравнения существования самолета определить геометрические параметры крыла, взлетную массу, массу конструкции, потребные характеристики силовой и энергоустановки. С помощью данной методики спроектирована эскизная модель БПЛА Прото-тип-2Э, спроектирована рабочая модель и изготовлен БПЛА Прототип-2Т. Проведено сравнение результатов расчетов аэродинамических характеристик БПЛА с результатами расчетов в пакетах OpenFOAM и ANSYS Fluent. Наблюдается хорошее соответствие результатов расчетов. Расчеты проводились на ресурсах суперкомпьютера Томского государственного университета СКИФ Cyberia Приведенная методика расчета первого приближения БПЛА может быть использована для расчета и экстраполяции кривых аэродинамических характеристик, летно-технических характеристик на область рабочего диапазона углов атаки и скоростей с поправками по результатам численного моделирования одной, двух или более точек. Также метод может быть использован для расчета начального приближения в задачах оптимального поиска с целью сокращения времени расчета.
Ключевые слова
беспилотный летательный аппарат,
проектирование БПЛА,
расчет первого приближения,
аэродинамические критерии эффективностиАвторы
Исмаилов Куат Кайратович | Томский государственный университет | младший научный сотрудник суперкомпьютерного центра | kuat@ftf.tsu.ru; mendikjan@gmail.com |
Всего: 1
Ссылки
Федоров Л.П., Михайлов Ю.С. Определение оптимальных режимов крейсерского полета высотного беспилотного летательного аппарата // Научный вестник Московского государственного технического университета гражданской авиации. 2013. № 2 (188)., С. 72-76.
Егер С.М., Мишин В.Ф., Лисейцев Н.К., Бадягин А.А., Ротин В.Е., Склянский Ф.И., Кон драшов Н.А., Киселев В.А., Фомин Н.А. Проектирование самолетов : учебник для вузов. М. : Машиностроение, 1983. 616 с.
Мизес Р. Теория полета / пер. с англ. А.Н. Рубашова. М. : Иностр. лит., 1949. 696 с.
Торенбик Э. Проектирование дозвуковых самолетов / пер. с англ. Е.П. Голубков. М. : Машиностроение, 1983. 648 с.
Пейгин С.В., Орлов С.А. Оптимальное аэродинамическое проектирование конфигурации «крыло-фюзеляж» широкофюзеляжного дальнемагистрального самолета // Вестник Томского государственного университета. Математика и механика. 2020. № 63. С. 115124. doi: 10.17223/19988621/63/10
Степанов К.А., Тимченко С.В. Аэродинамическое проектирование изолированного трехмерного крыла беспилотного летательного аппарата // Вестник Томского государственного университета. Математика и механика. 2018. № 54. С. 118-130. doi: 10.17223/19988621/54/10
Лойцянский Л.Г. Механика жидкости и газа. М., 2012. 210 с.
Шлихтинг Г. Теория пограничного слоя / пер. с нем. Г.А. Вольперт с 5-го нем. изд., исправленный по 6-му (амер.) изд., под ред. Л.Г. Лойцянского. М. : Наука, Гл. ред. физ.-мат. лит., 1974. 711 с.
ГОСТ 22833-77 Характеристики самолета геометрические. Термины, определения и буквенные обозначения. М. : Изд-во стандартов, 1987. 24 с.
Исмаилов К.К., Кагенов А.М., Костюшин К.В., Орлов С.А. Разработка цифровой модели БПЛА самолетного типа // 19-я Международная конференция «Авиация и космонавтика», 23-27 ноября 2020 г., Москва : тез. М. : ПеРо, 2020. С. 59-60.
Wilcox D.C. Turbulence Modeling for CFD. La Canada, CA : DCW Industries, 1993. 460 р.
Menter F.R. Zonal two-equation k-ю turbulence model for aerodynamic flows // AIAA. 1993. V. 93. Art. 2906. doi: 10.2514/6.1993-2906
Kagenov A.M., Kostyushin K.V., Ismailov K.K., Kostyushina N.O., Orlov S.A., Prokhanov S.A. The development of a cloud system for investigation of UAVs aerodynamic characteristics //j. Phys.: Conf. Ser. 2020. V. 1488. P. 1-5. doi: 10.1088/1742-6596/1488/1/012017
Патанкар С. Численные методы решения задач теплообмена и динамики жидкости : пер. с англ. М. : Энергоатомиздат, 1984. 152 с.