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

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

Том X 197 9 № 1

УДК 532.582.33

МЕТОД РАСЧЕТА ОСЕСИММЕТРИЧНОГО ОБТЕКАНИЯ ИДЕАЛЬНОЙ ЖИДКОСТЬЮ КОЛЬЦЕВОГО КРЫЛА С ЦЕНТРАЛЬНЫМ ТЕЛОМ

3. Г. Левшина, Л. А. Маслов

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

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

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

1. Постановка задачи. Течение идеальной несжимаемой жидкости вокруг системы контуров, образующих осесимметричные тела, моделируется движением системы источников и стоков, непрерывно распределенных по поверхности тела ш, и кольцевых вихрей, располагающихся внутри тела. Потенциал такого течения Ф(я, у, г, () существует везде вне замкнутых областей, окруженных контурами, и удовлетворяет уравнению Лапласа. Условие

непроницаемости поверхности контуров запишем в традиционной форме

<ЭФ

0-1)

где потенциал Ф включает в себя потенциалы слоя источников

—У —V

и вихрей, V—переносная скорость точек поверхности ш, п — внешняя нормаль к поверхности каждого контура.

С острых кромок областей, образованных контурами, в соответствии с гипотезой Чаплыгина — Жуковского в поток сходят

струи жидкости. Обозначив единичным вектором <7 направление, перпендикулярное линии схода струй, условие Чаплыгина — Жуковского можно использовать в виде равенства нулю компонента

—У

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

(V + V + йГ) ? = 0, (1.2)

где V — вектор скорости жидкости, индуцируемой слоем источников и стоков, т— вектор скорости, индуцируемой вихрями.

Для определения интенсивности СЛОЯ ИСТОЧНИКОВ (А условие (1.1) можно записать в виде интегрального уравнения

“Г*

2^(р) + ®я,

(1.3)

а условие (1.2) на задней кромке кольцевого крыла конфигурации, показанной на фиг. 1, в виде уравнения

IX к - — Ук д

Як д

(1.4)

здесь Р — произвольная расчетная точка поверхности ш, <3—переменная точка этой поверхности, = причем /?к — (}РК, где точка Рк принадлежит задней кромке крыла, Г — параметр, характеризующий циркуляцию вихревой системы крыла, 1адк — компонент скорости, индуцируемой вихревыми особенностями при Г = 1, вдоль —>■

вектора <7 в точке задней кромки крыла, \/к — вектор переносной скорости в точке задней кромки.

Фиг. 1

Совместным решением уравнений (1.3) и (1.4) определяется распределение интенсивности источников ^ и параметра циркуляции Г, после чего могут быть подсчитаны скорости жидкости в любой точке пространства.

Задача, описываемая уравнениями (1.3) и (1.4), является нелинейной, так как направление вектора определяется самим решением. Она может быть замкнута с помощью метода последовательных приближений, заключающегося в последовательном определении действительного направления линии схода струй с задней кромки профиля. При этом необходимо многократное решение уравнений (1.3) и (1.4), требующее значительного объема вычислений. В то же время изменение этого направления ограничено малым углом заострения задней кромки (см., например, [6]), а биссектриса этого угла, как правило, почти параллельна оси симметрии крыла.

В данной статье в целях сокращения времени решения задачи продольного обтекания кольцевого крыла линия схода струй с острой задней кромки задается параллельной продольной оси. При заданном направлении линии схода струй задача оказывается линейной относительно потенциала возмущенных скоростей. В данном —►

случае вектор является единичным вектором радиальной координаты цилиндрической системы.

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

2. Расчетные .формулы. Для записи расчетных формул используются результаты работы [5], представленные в несколько ином виде, удобном для использования в случае нескольких контуров. Для простоты рассматривается конфигурация, образованная тремя контурами (фиг. 1). Кольцевое крыло считается образованным двумя контурами. Параметрам, характеризующим внешний контур профиля, присваивается индекс 1, для внутреннего контура используется индекс 2. Центральное тело образуется одним контуром, его параметрам присваивается^ индекс 3. Каждый контур в цилиндрической системе координат, показанной на фиг. 1, представляется в параметрическом виде

лг = л:(5), г = г(«), (2.1)

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

х' = + дх/дя, г' = + дг/дз, (2.2)

где знак „ + “ используется в том случае, когда область течения при движении по контуру в сторону возрастания « остается слева, и знак „—когда справа.

Расчетной точке пространства присваиваются координаты д:, г, 0, а если точка принадлежит поверхности тела, то и значение параметра в. Аналогичные координаты текущей точки поверхности, по которой ведется интегрирование, обозначаются ?, р, ^ и в с индексами, соответствующими данному контуру. Продольные координаты вихревых колец и их радиусы также обозначаются буквами % и р.

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

здесь ^(о,) — приведенная интенсивность источников на соответ* ствующем контуре (V = 1, 2, 3), а ядра интегралов обозначены

К {к2) и Е{№) — полные эллиптические интегралы первого и второго рода с модулем £2 = 4гр„/7?|.

Для точек на оси симметрии (г = 0) интегралы (2.3) упрощаются аналогично [5].

Кольцевые вихри внутри крыла располагаются в соответствии с распределением расчетных точек на профиле крыла. При условии, что координаты х! расчетных точек на контурах внешней и внутренней поверхностей совпадают, эти кольцевые вихри имеют продольные координаты %} = Х] и радиусы, равные радиальным координатам точек срединной поверхности профиля ру = (гх} + г2;-)/2. Здесь, как и ранее, индексы 1 и 2 относятся к внешнему и внутреннему контурам, а номер расчетной точки /' меняется от 1 (носовая точка) до N (хвостовая точка). При этом в самих концевых точках вихри не ставятся, и для лучшего моделирования обтекания закругленного носка первый вихрь размещается в той точке с номером т, где половина местной толщины профиля равна или меньше отстояния точки от носка, исключая т= 1. Распределение циркуляции этих дискретных вихрей по хорде задается из условия, что они аппроксимируют вихревой слой, интенсивность которого пропорциональна местной толщине профиля с^ = (г1;- — г2;)/Ь, где Ь — хорда профиля.

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

з ^

4 = 1 О

з

(2.3)

(2.4)

А = -±~, б, = £(**), б,-К (**)-£(*»),

2

Я? = (х - ^)2 + (г - Р,)2, = (х - ЪУ + (г + Рч)2;

(2.5)

ввести относительные координаты £0у. и р0;- и безразмерный параметр циркуляции Г по формулам

Ьо) = (Х — Ро 1 = Г1?]> Г; = 4^; С, Д?у 1/оо Г, (2.6)

где Д^ = 5/+1 — 5у_1.

Опуская промежуточные выкладки, заимствованные из [8], для вычисления безразмерных компонентов скоростей, индуцируемых вихрями, получим формулы

О*, г) = Г®, , (X, г), (х, г) = Г®,, (X, г), (2.7)

где

N—1 N—1

X = X С1Ах] Кх)> г = 2 Дл7 КФ (2.8)

/=т /=т

(2.9)

здесь К. и £ — полные эллиптические интегралы первого и второго рода с модулем &2 = 4р0у-/й|.

Для точек на оси симметрии (г = 0) вместо формул (2.9) следует пользоваться выражениями

= К,1 = 0. (2.10)

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

их=\ + Ъх + ™х, иг = ъг+тг. (2.11)

Интегральное уравнение (1.3) после некоторых преобразований для каждого контура можно записать в виде

я = гг' + г [г' (да, + тх) - х' {уг + Юг)], (2.12)

где неизвестная интенсивность источников g входит под знак интеграла благодаря формулам (2.3), а неизвестный параметр циркуляции Г оказывается в правой части благодаря формулам (2.7). Если воспользоваться компонентами относительной скорости (2.11), то уравнение (2.12) получит более компактный вид

ё = г{г'их — х'иг). (2.13)

Условие (1.4) для определения параметра циркуляции получается в виде

TwtrK = — vrЮ (2.14)

где единичная скорость от вихревой системы на задней кромке определяется формулами (2.8), а неизвестная интенсивность источников g оказывается в правой части под интегралом в соответствии с (2.3).

Относительная скорость на поверхности тела представляется в виде касательного компонента, равного

их = х' их + г' иг, (2.15)

■а безразмерный коэффициент давлений вычисляется с помощью •формулы Бернулли р = 1—и2г.

3. Методика вычислений. Решение интегрального уравнения (2.13) по аналогии с [5] определяется в конечном числе дискретных точек каждого контура, называемых расчетными.

Вычисление интегралов (2.3) в случае нескольких контуров выполняется в зависимости от расположения расчетной точки относительно контура. Если расчетная точка отстоит от поверхности, образованной контуром с номером v на расстоянии г0, превышающем 5% длины дуги Lv этого контура, то интегрирование по нему ведется при помощи обычного правила трапеций по заданным точкам контура. В противном случае интегралы по такому контуру вычисляются как несобственные при помощи замены переменных

0,5L ЛЗ

О - а0 = sign (/г) _ h2................., (3.1)

где а0 — координата s той точки контура, расстояние которой от расчетной точки минимально.

При этом замена (3.1) используется лишь на участке контура

I — з0 I < 0,05 U (3.2)

аналогично тому, как это делается в работе [5] для одного контура.

Интегральное уравнение (2.13) вместе с (2.14) решается методом последовательных приближений. Полагая Г = 0 и в качестве начального приближения g-(0) = 1/2 гг', по формулам (2.3) подсчитывается vj,°l в точке задней кромки крыла, а затем по формуле (2.14) значение параметра циркуляции Г<°> в нулевом приближении. В первой итерации сначала вычисляются в расчетных точках всех контуров значения интенсивности g**1) по формуле (2.13), где необходимые соотношения (2.11) используют g(0> и Г<°). Интенсивность в первом приближении полагается равной

£(,) = 4-(£(0) + £*(1))- (3-3>

Затем вместо уравнения (2.14) решается уравнение

дг«'1,к = -и,к, (3.4)

где при вычислении игК с помощью (2.11) используются и Г<°>. Параметр циркуляции в первом приближении определяется равенством Г<!) = Г(°) + ДГО). Для итерации с номером п последовательность вычислений выглядит следующим образом: по формуле (2.13) определяется g* <") по значениям и Г*"-1), подсчитывается

g-(n> = l/2(g(n-1) -f- g* <">), по формуле (3.4) подсчитывается прираще-

ние параметра циркуляции ДРл>, используя gW и Г*"-1), после чего определяется новое значение параметра циркуляции:

г<«) = Г(л-!) + ДГ<л>.

На этом этапе выполняется оценка сходимости процесса итераций-Достаточно сравнить между собой лишь интенсивность источни' ков и стоков и продолжать итерации до тех пор, пока

] g(n)_g(n-l) |max>s I rr' |max,

где считалось достаточным положить е = 0,005.

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

Информация о рассчитываемой конфигурации задается набором таблиц координат х и г расчетных точек для каждого контура.

4. Результаты расчетов. Изложенная методика реализована в программах на языке ФОРТРАН. В случае циркуляционного обтекания изолированного кольцевого крыла или с центральным телом не представляется возможным сравнение с точными решениями. Для проверки правильности работы программ вычислялось распределение скоростей по поверхности кольцевого крыла с профилем Жуковского, когда его диаметр составляет десять хорд и обтекание профиля можно считать плоским. Результаты этого расчета удовлетворительно совпали с точным решением для этого профиля под нулевым углом атаки в плоском потоке.

На фиг. 2 в качестве иллюстрации построены результаты расчета безразмерного коэффициента давлений р на внешней (кривые 1) и внутренней (кривые 2) поверхностях кольцевого крыла с симметричным профилем Жуковского толщиной 12% и на поверхности центрального тела вращения (кривые 3), образованного тем же профилем. Диаметр кольцевого крыла равен хорде и длине центрального тела. Сплошные линии соответствуют давлениям на кольцевом крыле с центральным телом. Для сравнения пунктирными линиями построены давления для изолированного кольцевого крыла (кривые / и 2), а также центрального тела в безграничном потоке (пунктирная кривая 3). Можно отметить, что влияние центрального тела в этом примере почти не сказалось на обтекании внешней поверхности и привело к некоторому увеличению разрежения на внутренней поверхности крыла. Заметно изменилось разрежение на центральном теле по сравнению с обтеканием безграничным потоком.

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

Фиг. 2

р

-0,3

-0,2

-0,1

О

0,1

0,2

0,3

Фиг. 3

■Р

-0,8

-0,1

-0,2 О 0,2 о

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

На фиг. 3 и 4 выполнено сравнение расчетных давлений с экспериментальными для двух кольцевых крыльев без центрального тела, полученных в работах соответственно [9] и [10]. Сплошными линиями построены результаты расчетов, точками—-экспериментальные значения коэффициента давлений. Светлые точки соответствуют внешнему контуру, закрашенные — внутреннему. Профили крыльев построены на этих же фигурах. Первое крыло имеет симметричный профиль ИАСА-бб толщиной 6% и диаметр срединной поверхности, равный хорде. Профиль второго крыла (см. фиг. 4) имеет сложную форму, полученную в работе [10] в результате

Фиг. 4

вычислений по заданному распределению особенностей, толщину 6%, а наименьший внутренний диаметр крыла равен 1,27 хорды.

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

ЛИТЕРАТУРА

1. Сидоров О. П. Обтекание простейшего кольцевого крыла. Труды КАИ, вып. 29, 1955.

2. Trulin D. J. A method for computing the pressure distribution about an axisymmetric nacelle In subsonic flow. Ames, 1968 (Dissertation).

3. Hess J. L., Smith A. M. O. Calculation of potential flow about arbitrary bodies. .Progress in Aeronautical Sciences", vol. 8, 1967.

4. Y о u n g C. A computer program to calculate the pressure distribution on an annular aerofoil. ARC CP N 1217, 1972.

5. Маслов Л. А. Метод расчета обтекания тела вращения любой формы при произвольном движении в идеальной жидкости. „Ученые записки ЦАГИ“, т. 1, № 2, 1970.

6. Mangier К. W., Smith J. Н. В. Behaviour of the vortex sheet at the trailing edge of a lifting wing. , Aeron. J.“, XI, vol. 74, N719. 1970.

7. Hess J. L. Review of integral-equation techniques for solving potential-flow problems with emphasis on the surface-source method. „Comput. Meth. Appl. Mech. and Engin.“, vol. 5, N 2, 1975.

8. Белоцерковский С. М. Тонкая несущая поверхность в дозвуковом потоке газа. М. „Наука*, 1965.

9. Hough G. R. The aerodynamic loading on streamlined ducted bodies. „Developments in Mechanics. Conference held at Case Inst, of Technology, April, 1963“, vol. 2, part 1, 1965.

10. R у a 11 D. L., Collins 1. F. Design and test of a series of annular aerofoils. .ARC R and M“, N 3492, 1967.

Рукопись поступила 20(11 1978 г.