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

УЧЕНЫЕ ЗАПИСКИ Ц А Г И Т о м VII 197 6

М I

УДК 533.6.011.32

К РАСЧЕТУ ОБТЕКАНИЯ ФЮЗЕЛЯЖА ПРОИЗВОЛЬНОЙ ФОРМЫ ПРИ МАЛЫХ СКОРОСТЯХ

Л. А. Маслов, В. П. Юшин

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

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

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

В работе[3] предложен метод расчета обтекания трехмерного тела, отличающийся от метода [1] способом решения основного интегрального уравнения. Не уступая в точности методу [1], метод

[3] оказывается менее трудоемким и позволяет использовать для расчетов машины среднего быстродействия (до 1 часа на машине типа М-20).

В данной работе совершенствуется методика [3]. Вместо цилиндрических координат, которые несколько ограничивали допустимую

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

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

2*1* (PH ■ (1.1)

О)

здесь Р — произвольная расчетная точка; <3 — текущая точка поверхности со; /? = <ЗР; п — нормаль к поверхности ш, направленная

внутрь жидкости; ^(Р) —скорость точек поверхности тела (переносная скорость).

Уравнение (1.1) является интегральным уравнением Фредгольма второго рода, имеющим единственное решение в классе непрерывных функций, если поверхность <0 принадлежит к классу Ляпунова

[4]. Ядро уравнения имеет слабую особенность: < С/7?2-а при

Я -*■ 0. Здесь С — постоянная, <х>0— показатель Ляпунова поверхности со. Особый интеграл уравнения существует.

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

После определения из (1.1) интенсивности ^ скорость жидкости относительно фюзеляжа в любой точке его поверхности будет найдена по формуле

£{Р) = 2«рп + ^г(а)-^<1и-£(Р). (1.2)

ш

Если точка Р не принадлежит о>, в правой части (1.2) следует опустить первое слагаемое.

2. Представление поверхности. Используется связанная с телом

* ►—»-

система декартовых координат хуг с ортами гу'£. Координатная плоскость хОу является вертикальной плоскостью симметрии (фиг. 1). Продольная ось х направлена от носа фюзеляжа к его хвосту. Положение оси х по высоте фюзеляжа произвольно. Координатная плоскость уОг должна проходить через носовую точку фюзеляжа. Расстояние вдоль оси х от начала О до точки пересечения х с плоскостью, параллельной уОг и проходящей через хвостовую точку, считается длиной фюзеляжа Ь. Обводы фюзеляжа задаются набором поперечных сечений, перпендикулярных оси х.

При вычислениях в каждом поперечном сечении х — const на поверхности фюзеляжа вводится координата длины дуги s, отсчитываемая от верхней точки, лежащей в плоскости симметрии хОу. Направление отсчета s в окрестности этого начала совпадает с направлением боковой оси z. Значение параметра s в нижней точке поперечного сечения в плоскости хОу, т. е. полупериметр сечения,

обозначается буквой I. Эта величина используется для представления местной координаты 5 в безразмерном виде «= 5//, диапазон изменения которой для всего сечения оказывается в пределах от О до 2.

Обводы поперечного сечения представляются в параметрическом виде

у=у{8), х = г(з). (2.1)

Для выражения дифференциальных элементов поверхности введены обозначения

ду

Pi

ds

дг

Р2 = ИГ

Рг-

дг

ду

(2.2)

При вычислении рх и рг (я = const) размерности s, у и z приняты одинаковыми. Это дает удобное соотношение р\-\-р\=\- Производные ду/дх и дг/дх вычисляются при постоянных значениях безразмерной координаты s.

В соответствии с формулами дифференциальной геометрии элемент поверхности представляется в виде

du> = I (х) Y1 ■+■ р\ dx ds. (2.3)

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

я = (1 +pl)~v2(p3i+pj — Pi к); т = (1 +р1)~Щ ,(i—'PiPd+PiPtb)-, (2-4)

b=PU + Ptk,

где п — нормаль, вектор Ь направлен по линии пересечения касательной и поперечной плоскостей, вектор х образует с п и Ь правую систему координат и имеет положительную проекцию на ось х.

Форма поперечных сечений и их распределение относительно оси х в данном случае совершенно произвольны. Как и при использовании метода [3], не должны иметь места касательные к поверхности плоскости, перпендикулярные оси х (рв ограничено). Исключением являются концевые точки фюзеляжа, где указанные плоскости должны иметь единственную точку касания.

Форма фюзеляжей современных дозвуковых самолетов удовлетворяет этому требованию. В концевых точках, где производная Ра~>-оо и полупериметр 1(х)^ 0, необходимо вычислить предел

с («) == Нт [I (х) рь (х, 5)]. (2.5)

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

у (х, 7) = ах (1) х+д1 (I) Ух + ук, , ^ ^

z (х, s) = аг (s) х + b2 (s) Yх.

Тогда х

с = 4-(Mi-Mi). = = (2.7)

здесь хк и ук — координаты концевых точек.

Для хвостовой точки в (2.6) необходимо х заменить на L — х. Определение коэффициентов Ьл и Ь2 и их производных осуществляется по значениям координат двух расчетных сечений х = const, ближайших к концевой точке.

3. Расчетные формулы. Расчетной точке пространства Р присваиваются координаты х, у, z, а если эта точка лежит на поверхности фюзеляжа, то значение параметра s. Координаты текущей точки Q поверхности о> обозначаются соответственно 5, тг), С и а. Вектор

#=(*-6)7+ (у - ч) У + (* — С) Л, (3.1)

а его модуль

/? = /(*- + (jf - т,)2 + {z- С)2. (3.2)

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

v = — {i cos а + у sin a) Voo, (3:3)

где Voo—скорость связанной с телом системы координат относительно некоторой неподвижной системы.

Вместо интенсивности у- удобно пользоваться функцией g, связанной С f* соотношением

£ = 2*/УТ + Лг-- (3-4)

vco

Проекции безразмерной возмущенной скорости жидкости запишем в виде

1 2

Vi = J j SKi da d\, i = x, у, z, (3.5)

0 0

где

2л/?з ’ У 2л/?3 ’ 2kR? ■

Интегральное уравнение (1.1) для определения функции g при помощи (2.3), (2.4) и (3.1) — (3.5) преобразуется к следующему виду: g = I (—/?з cos а - р2 Sin а — р3 vx — p2vy + рг v2), (3.6)

где неизвестная g в правой части содержится в выражениях для vx, vy и vz (3.5).

Аналогичным образом получаются формулы для проекций безразмерных касательных скоростей на поверхности фюзеляжа Ит = (1 +Pl)-m(cosa-pip3sina + vx-p2psvy-\-p1pBvz)- 1

«г- = Pi sin« +А»у + Рг^г I

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

их = cos а -|- vx, иу = sin а-f- vy, uz — vz. (3.8)

Для концевых точек интегральное уравнение (3.6) с учетом (2.5) преобразуется к виду

i") = c( ”, s) - cosa —, (3.9)

а при определении относительной скорости достаточно вычислить компонент «у. Давления подсчитываются при помощи формулы Бернулли.

4. Методы вычислений. Все вычислительные приемы здесь оставлены теми же, что и в работе [3]. Определение интенсивности источников при решении уравнения (3.6) и последующий подсчет скоростей и давлений осуществляются в конечном числе дискретных расчетных точек поверхности фюзеляжа. Последние располагаются в произвольно выбранных расчетных сечениях х-= const. В программе для машины БЭСМ-6 может использоваться до 63 таких сечений. На одной половине контура сечения располагаются произвольным образом до 50 расчетных точек. Общее число этих точек на всем фюзеляже не должно превышать 1280. Информация о расчетных точках, являющаяся информацией о поверхности фюзеляжа, вводится в ЭЦВМ в виде таблиц ординат. Увеличение числа расчетных точек значительно облегчает рассмотрение фюзеляжей довольно сложной формы.

Интегральные уравнения решаются методом последовательных приближений. Решение считается найденным, если

I ёп. Sn—1 Imax _ /л i\

<е, (4.1)

| I (р3 COS а-\-р2 sin

a) I max

где n — номер приближения и г = 0,001.

Несобственные интегралы (3.5) вычисляются при помощи замен переменных

(4-2)

(4-3)

Эти замены соответствуют заменам, приведенным в работе [3], с той лишь разницей, что последняя формула приспособлена к интегрированию по контуру а. Параметр о0 характеризует точку на контуре ? = const, расстояние которой г от расчетной точки,

измеренное в поперечной плоскости, минимально. Замена (4.3) используется только при условиях

|*-е|<0,05, | а~ — о~01 -< 1/12, г<0,064. (4.4)

Вне этих пределов внутренний интеграл формул (3.5) вычисляется непосредственно по переменной о по правилу трапеций. Для интегрирования по переменным Лиф используется правило трапеций

Фиг. 2

с постоянным шагом, который выбран тем же, что и в работе [3]: ДЛ = 384~1/2 и Дф = 2-4. Необходимые в процессе такого интегрирования промежуточные значения координат и интенсивности источников между заданными расчетными точками определяются двойной квадратичной интерполяцией. _

Ядра интегралов (3.5) при % = х, а также на линии о = о0 при условиях (4.4) из сумм квадратурной формулы исключаются, кроме концевых точек, когда они являются расчетными. В этом случае в сумме для скорости <ох учитывается предельное значение, полученное при помощи (2.6) и (4.2) и равное для носовой точки

К "г") = ГШ(Ь21 + Ь^2 ’

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

5. Пример расчета. Для проверки точности изложенной методики на машине БЭСМ-6 выполнены расчеты скоростей на поверхности трехосного эллипсоида с соотношениями полуосей вдоль х, у иг соответственно 4:2:1. Расчетные значения скоростей при продольном обтекании показаны точками на фиг. 2 для вертикальной плоскости симметрии г = 0 (светлые точки) и для горизонтальной плоскости симметрии _у = 0 (черные точки). Сплошные линии соответствуют известным точным значениям скоростей на поверхности этого эллипсоида. Получено хорошее соответствие точных и расчетных значений.

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

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

причем для удобства построения эпюры давлений поперечные размеры увеличены по сравнению с продольными в 2,5 раза. Значения р отложены по нормали к контуру: разрежение — по внешней нормали, давление — по внутренней. Сплошными линиями показаны расчетные значения (при числе М = 0), точками — результаты эксперимента, полученные в аэродинамической трубе при числе Ке = 2,1Х107 и М = 0,15 (закрашенные точкй) и М = 0,4-(светлые точки).

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

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

1. Hess I. L., Smith A. M. O. Calculation of potential flow about arbitrary body shapes. International Symposium on Analogue and Digital Techniques Applied to Aeronautics, Li6ge, 1963, Bruxelles, 1964.

2. Хоиичев В. И. Численное решение задачи о потенциальном обтекании пространственных тел. Изв. Сиб. отд. АН СССР, серия техн. наук, вып. 2, 1973.

3. Маслов Л. А. Произвольное движение продолговатого тела в идеальной жидкости. „Изв. АН СССР, МЖГ“, 1966, № 6.

4. М и х л и н С. Г. Лекции по линейным интегральным уравнениям. М., Физматгиз, 1959.

Рукопись поступила 311 1975