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

Том XXXVIII

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

№ 1 — 2

УДК 533.6.011.3/.5:532.582.2 533.6.011.3/.5:629.7.024.36

ОБТЕКАНИЕ ЭЛЛИПТИЧЕСКОГО ЦИЛИНДРА И ПЛАСТИНЫ

В ПРИБЛИЖЕНИИ ОЗЕЕНА

А. В. КАШЕВАРОВ

На основе точных решений уравнений Озеена, предназначенных для описания медленных движений вязкой жидкости при числах Рейнольдса Re << 1, рассчитаны картины стационарного симметричного обтекания эллиптического цилиндра и пластины вплоть до числа Re = 10. Найдены критические значения Яе*, соответствующие началу образования вихревой зоны за эллипсом, в зависимости от соотношения его осей и ориентации.

Определены параметры вихревой зоны при различных Re, а также получены значения коэффициента сопротивления цилиндра и пластины.

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

в то время картины обтекания этих тел при Яе ~ 1 на основе полученных сложных аналитических выражений в виде рядов произведений различных специальных функций. Было лишь выяснено, что решения уравнений Озеена дают качественно верные картины обтекания, предсказывая появление вихревой зоны за телом при некотором критическом числе Яе*. Для кругового цилиндра в [4] было найдено Яе* = 1.51. Позднее с помощью ЭВМ были рассчитаны линии тока вокруг сферы [5] и кругового цилиндра [6] в приближении Озеена на основе имеющихся общих аналитических решений, но постоянные интегрирования в них определялись приближенно в результате минимизации квадрата скорости и2 течения на

контуре С тела, т. е. искали шт еи 2ё1. В случае сферы было найдено критическое число

с

Яе* = 3.81.

С появлением мощных ЭВМ интерес к аналитическим решениям уравнений Озеена практически иссяк, и основное внимание стало уделяться численным решениям уравнений Навье —Стокса. С целью сравнения приближенного и точного подхода к описанию движения вязкой жидкости в [7] вновь вернулись к рассмотрению вопроса обтекания тел в приближении Озеена при Яе ~ 1, сосредоточившись на параметрах вихревой зоны за телом. В отличие от [5, 6], картины течения вокруг сферы и кругового цилиндра в [7] были рассчитаны на основе точных решений [1], где постоянные функциональных рядов определяются из условий торможения потока

на поверхности тела. Были подтверждены результаты предыдущих исследований относительно критического числа Яе*. Сравнение показало, что вихревая зона образуется при более низких числах Яе*, чем следует из уравнений Навье — Стокса, и является более протяженной и широкой при одном и том же числе Яе > Яе*. Особенно большие расхождения имеются для пространственного обтекания сферы, однако в случае плоского обтекания кругового цилиндра

оказалось, что приближение Озеена достаточно хорошо предсказывает ширину вихревой зоны. Так, при Яе = 15 угловая координата отрыва лишь на 6% отличается от значения, полученного из численного решения уравнений Навье — Стокса.

Можно предположить, что согласие между приближенным и точным описанием вязкого течения улучшиться для вытянутого эллиптического цилиндра, обтекаемого вдоль большой оси, так что в этом случае точные решения приближенных уравнений Озеена могут конкурировать с приближенными численными решениями точных уравнений движения вязкой жидкости. Нужно сказать, что расчет обтекания эллиптического цилиндра в рамках полных уравнений Навье —Стокса, видимо, представляет собой более трудную задачу, чем в случае кругового цилиндра, так как число работ, посвященных этой проблеме, значительно меньше. Стационарное течение вокруг эллиптического цилиндра, обтекаемого вдоль большой оси, рассчитано в [8]. Расчет обтекания тонкого эллипса под углом атаки проведен в [9]. Однако детальное численное исследование

параметров вихревой зоны за эллиптическим цилиндром отсутствует. Исключение составляет работа [10], в которой рассматривалось обтекание плоской пластины, установленной перпендикулярно набегающему потоку. Численное решение уравнений Навье — Стокса для продольного обтекания пластины выполнено в [8, 11].

Продольное и поперечное обтекание пластины в приближении Озеена рассматривалось в [12] и [13, 14] соответственно. Рассмотрение проводилось на основе решения в виде интегрального уравнения для распределения фундаментальных особенностей на поверхности пластины, которое решалось приближенно. В [12] внимание было сосредоточено на поведении коэффициента трения вблизи задней кромки. Зависимость коэффициента сопротивления поперечно обтекаемой пластины от числа Яе получена в [13]. Лишь в [14] представлено описание вихревой зоны за пластиной в приближении Озеена.

Насколько известно автору, до сих пор не было получено описание вихревой зоны за эллиптическим цилиндром на основе точных решений уравнений Озеена. Этот вопрос является предметом рассмотрения данной статьи. Рассмотрение ограничено только случаем симметричного обтекания [2].

1. Основные уравнения и формулы. Хорошо известные линеаризованные уравнения Озеена движения вязкой несжимаемой жидкости могут быть записаны в векторной безразмерной форме следующим образом:

ди 1

— = -V» + Ди,

дх 2к

и = 0.

(1.1)

(1.2)

Здесь и — поле возмущений скорости, обезразмеренное через скорость набегающего потока Ц», направленную вдоль оси х декартовой системы координат (х, у); р — возмущение давления,

отнесенное к рЦ» (р — плотность жидкости); к = Яе/2. В качестве характерной длины, определяющей число Яе, возьмем полусумму большой и малой полуосей эллипса (а + Ь)/2, так что

в пределе кругового цилиндра характерным размером будет его радиус, а в другом предельном случае плоской пластины — четверть ее длины Ы4.

Уравнения (1.1), (1.2) удовлетворяются при

1 д и = ^Ф +—Vx-xi, р = —.

Здесь 1 = (1, 0) — единичный вектор, а функции х и ф в свою очередь удовлетворяют уравнениям

ДХ-2кдХ = 0, (1.3)

дх

Вихрь ю равен

ю = rot u = -V% х i. (1.4)

Решение задачи проводится в эллиптической системе координат (S, п), связь которых с обезразмеренными декартовыми в случае обтекания вдоль большой оси эллипса дается соотношениями:

x = 2exp(-^0)ch ^cos п, y = 2exp(-^0)sh ^sin п.

Здесь ^0 — радиальная эллиптическая координата контура эллипса. Предельные случаи кругового цилиндра и пластины соответствуют ^0 ^ да и ^0 = 0.

Общее решение уравнения (1.3) в эллиптических координатах имеет вид:

да

X = exP (2Kch £, cos п) X Pm Fekm (S -q)cem (П - q)-

m=0

Здесь cem(n, -q) — четные периодические решения обычного уравнения Матье; Fekm(S, -q) — непериодические вторые решения модифицированного уравнения Матье, стремящиеся к нулю при S ^ да; q = К2 = k2exp(-2^0); pm — постоянные интегрирования, подлежащие определению из условия равенства нулю скорости на поверхности цилиндра.

Повторяя ход решения [2], но войьшей степени используя дости жения теории функций Матье [15], придем к следующему выражению для функции тока (штрихом обозначена производная по S):

у = 2sh (S - S0 ) sin п +

" sin пп

££|"Fekm (^ - q) am, n (S)- e n(S S0) Fekm (^ - q) am, n (S0 )

n=1 m=0

e S0 да да , г (15)

' —Z Z Pm {shS Fekm (^ - q)[am,n+1 (S) + am,n-1 (S)] -n=1 m=0

-e n(S S0) shS0Fekm (^ - q)\_am,n+1 (^0 )+ am,n-1 (S0 )]}

sin пп

n

Первый член правой части описывает функцию тока потенциального течения, остальные учитывают влияние вязкости.

В (1.5) функции ат, „(Е) представляют собой коэффициенты разложения функции х в ряд Фурье. Их можно найти, используя интегральное уравнение

П

Сет (Е - Ч) = У т | ехР (2ксЬ Е с08 П) сет (п - Ч) ^П- (1-6)

0

Здесь Сет(Е, -4) — первые решения модифицированного уравнения Матье; ут — собственные значения, которые равны

V2т = се2т (0, Ч)/пА0^ , V2т+1 = Эе2т+1 (0 Ч)/ ПкВ1 ^ ,

где а)2т), В 2т+1) — коэффициенты рядов представлений функций Матье се2т(п, ч), se2m+1(n, ч) по тригонометрическим функциям

ce2m (n q) = X A22m) cos 2Г П

r=0

se2m+1 (n, q) = Z 522 m1+1) sin (2r+1)п.

r=0

Из (1.6) сразу имеем

am,0 = 2Cem (S - q)/nVm .

Коэффициенты п-го порядка могут быть получены n-кратным дифференцированием (1.6) по S. Из-за громоздкости их выражения здесь не приводятся. Полезны рекуррентные соотношения:

am, 0 am,п

am,1 = » ГТ, am, n+1 = jTT-am, n-1. (1.7)

2k sh S к sh S

Условия прилипания приводят к следующей системе для определения постоянных вт-

1да

ТГ X Pmam, n (Fekm + n Fekm ) -

2k m=0

да f-2 n = 1 (18)

- T XPm Fekm (1 + 3e“2^° ) am, n+1 - (3 + e“^° ) am, nV , П

4

m=0

0, n = 2, 3,...

Значения функций Бект(Е, -4), коэффициентов ат, „(Е) и их производных по Е берутся в точке Е = Е0-

Коэффициент сопротивления цилиндра равен

CD = - { d= R" X Pm (Fekm am, 0 - Fekm am, 0 )•

f=fo Re m=0

0

Отметим, что CD здесь определяется как

С = D

CD -

ри2 (а + Ь)/2’

где О — сила сопротивления.

В случае обтекания эллипса вдоль малой оси эллиптические координаты вводятся по формулам:

X = 2ехр(-Ео)8Ь^С08П, У = 2ехр(-^о)сЬ^8ІПП.

Функция тока дается выражением вида (1.5), в котором первый член не изменяется, а далее все гиперболические синусы заменяются косинусами, оставшиеся величины — величинами со звездочками, причем

^к2т = ^к2т , а ^к2т+1 = ®ек2т+1>

а* _ т 2А2т)СЄ2т (Е, - Я) * ( 1)т+12кА1(2т+1) 8Є2т+і (Е, - Я)

а2т,0 V 1 ( /о \ ’ а2т+1, 0 V 1 г ( \ ’

Се2т ( V2 Я) Се2т+1 (П/2 Я)

» ( a*m ,0)

m,1 2кch E ,

к ch E

Здесь Оек2т+1 — непериодические вторые решения модифицированного уравнения Матье, соответствующие функциям 8е2т+1.

Постоянные вт определяются из решения системы

2. Процедура вычислений и результаты. Расчет функций тока по (1.5) представляет собой трудоемкую и кропотливую задачу. Она была решена на персональной ЭВМ типа Репйиш II с помощью программы, написанной на алгоритмическом языке БОЯТЯЛК.

с использованием библиотечной программы расчета собственных значений. Для расчета модифицированных функций Бесселя, определяющих функции Матье Fem, Cem, применялись представления функций K0, K\ в виде рядов и рекуррентные соотношения. Чтобы не происходило потери точности, суммированием соответствующих рядов получали значения двух функций Бесселя In, In-\ высокого порядка, значения этих функций более низких порядков вычислялись рекуррентно.

Для получения аналитических выражений коэффициентов ат, n(E) использовался пакет MAPLE V. Как следует из соотношений (1.7), при вычислении ат, n(E) порядков n > 0 возникает неопределенность вида 0/0, если Е ^ 0. Поэтому при расчете продольного обтекания плоской пластины, когда Е0 = 0, громоздкие и не приводимые здесь аналитические выражения для ат, n(0) находились в результате раскрытия неопределенности по правилу Лопиталя.

При небольших Е вычисления ат, n(E) высоких порядков n происходят с потерей точности. Поэтому расчет течения вокруг эллиптического цилиндра, обтекаемого вдоль большой оси, пришлось проводить с квадратичной REAL*16 точностью, обеспечивающей 32 десятичных знака, что позволило рассчитать обтекание цилиндра при Е0 - 0.4. Для реализации таких вычислений

на 32-разрядной ЭВМ использовался пакет программ [16]. Собственные значения и коэффициенты функций Матье уточнялись по методу [17]. Решение системы (1.8) осуществлялось с помощью библиотечной программы.

Расчеты обтекания эллиптического цилиндра проведены при нескольких значениях Е0. Прежде всего интересовались критическим числом Re* начала образования вихревой зоны. Из (1.4) имеем выражение для распределения интенсивности ю вихря на поверхности цилиндра:

Это выражение оценивалось в окрестности точки п = 0 путем его разложения в ряд Тейлора. Критерием появления вихревой зоны была смена знака главного члена разложения с ростом числа Яе. Чтобы уловить смену знака, был необходим учет членов ряда высокого порядка малости. Так, при Яе = 5 учитывалось двенадцать, а при Яе = 10 — шестнадцать членов ряда.

4 Ё P*m Fek*m Г(l - 3e“2Eo )am, n+1 - (з - e“2Eo)am, n-1

m=0

-2, n = 1 0, n = 2, 3,...

Вычисления коэффициентов a2 2m), B2+l+1) функций Матье осуществлялись по методу [15]

ю = -2е E0 ёР" [Fekm (^ -q) cem (n -q)ch E0sin n +

m=0

+Fek m (^ - q) cem (n - q)sh E0cos n].

Рассчитанные значения Re* при различных Е0 приведены в табл. 1. В колонке а представлены числа Re* в случае обтекания вдоль большой оси, а в колонке Ь — вдоль малой. Последний случай более легок для расчета, так как вычисления можно производить с двойной точностью.

Таблица 1

Е0 3 2 1 0.8 0.6 0.5 0.4 0.3 0.2 0.1

Яе, а Ь 1.53 1.50 1.60 1.43 2.34 0.99 2.92 0.81 4.15 0.59 5.34 0.47 7.47 0.36 0.25 0.15 0.06

Видно, что при Ео = 3, когда эллиптический цилиндр близок к круговому, значения Re* близки к 1.51. С уменьшением Ео критическое число Re* увеличивается при продольном обтекании и уменьшается при поперечном. Продольное обтекание пластинки происходит без образования вихря. При поперечном ее обтекании, как отмечалось ранее [14], вихревая пара образуется при произвольно малом числе Яе.

Примеры обтекания эллиптического цилиндра вдоль большой оси эллипса при Ео = 0.6 представлены на рис. 1, а (Яе = 5) и б (Яе = 10). Линия тока, примыкающая сзади к поверхности цилиндра, определяет границу вихревой зоны (у = 0). Остальные линии тока следуют с шагом 0.1 по у. Аналогичные картины обтекания эллипса вдоль малой оси приведены на рис. 2, а, б. На рис. 3, а, б показаны линии тока для продольного обтекания пластины (Е0 = 0) при тех же числах Яе. Поперечное обтекание пластины представлено на рис. 4, а (Яе = 0.5) и б (Яе = 1).

На рис. 5 приведены зависимости длины вихревой зоны Д, от числа Яе для различных Е0, а на рис. 6 — угла отрыва п (обозначения введены на рис. 2). Штрихпунктирные линии относятся к случаю кругового цилиндра (Е0 ^ да). Цифрами без штрихов обозначены зависимости для обтекания вдоль большой оси, со штрихами - вдоль малой. Кривые 1, 1' — соответствуют Е0 = 1; 2, 2' — 0.8; 3, 3' — 0.6; 4, 4' — 0.4; 5' — 0.2; 6' — 0.

Точками на рис. 5 представлены результаты [10], полученные из численного решения уравнений Навье — Стокса для обтекания установленной перпендикулярно потоку пластины. Протяженность вихревой зоны за пластиной с ростом Яе увеличивается быстрее в приближении Озеена. Что касается угла отрыва, то как уравнения Озеена, так и Навье — Стокса дают п = 90° для пластины при любых Re.

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

Численное интегрирование уравнений Навье — Стокса в случае обтекания эллипса

а)

б)

а)

вдоль большой оси проводилось в [8] при различных соотношениях его осей Ь/а, и числе Рейнольдса ReЬ = 40, определенного по длине большой оси Ь = 2а. В табл. 2 сравниваются данные настоящей работы и [8], последние пересчитаны по следующему правилу:

Со = 4СЬь/(1 + Л^о), Ео = агеШ(Ь/а).

Здесь СоЬ — исходное значение коэффициента сопротивления в [8], определенного

по длине большой оси СоЬ = о/рЦ^Ь . Числа Рейнольдса Яе и ЯеЬ связаны соотношением

Яе = ЯеЬ(1 + ЛЕо)/4. (Для получения значения Со при Е0 = 0.05 потребовались расчеты с 64 десятичными знаками).

Т аблица 2

Ь/а 0 0.05 0.1 0.2 0.4

Ео 0 0.050 0.100 0.203 0.424

Яе 10 10.5 11 12 14

Рис. 4

4 6

Со, Озеен 1.611 1.667 1.756 1.934 2.266

Со, Н - С [18] 1.292 1.295 1.302 1.323 1.377

Сравнение продолжено на рис. 7, где представлены зависимости Со(Яе). Штрихпунктирная кривая соответствует озееновскому обтеканию кругового цилиндра (Е0 ^ да), кривые 1 — 3 описывают зависимости при обтекании вдоль большой оси эллиптических цилиндров с Е0 = 1, 0.5 и 0 (пластина) соответственно. Кривые, обозначенные цифрами со штрихами, относятся к обтеканию цилиндров и пластины вдоль малой оси. Светлыми кружками приведены данные [18], найденные численным интегрированием уравнений Навье — Стокса в случае кругового цилиндра, а темными кружками — результаты [8] для продольного обтекания пластины конечной длины. Штриховая кривая соответствует закону Блазиуса для сопротивления продольно обтекаемой

—1/2

полубесконечной пластины, который описывается известной формулой [19] Соь = 1.328Яе_оЬ .

2

_- - 1_________________

-2-1 0 I 2x3

2 4 6 8 Яе 10

10

\W

\ ^ \ 0 4 ■ ° " — ^ O • 2 о 3 •

Рис. 8

Рис. 7

Сравнение показывает, что в приближении Озеена коэффициент С0 оказывается выше, чем это следует из расчетов вязкого обтекания в строгой постановке. При фиксированном Ео * 0 различие между значениями С0, полученными двумя способами, увеличивается с ростом Яе, достигая 73% при Яе = 10 для кругового цилиндра. С уменьшением Е0 согласие улучшается, и для пластины различие составляет всего лишь от 12% при Яе = 2.5 до 25% при Яе = 10. С ростом Яе коэффициент С0 пластины приближается к блазиусовской кривой.

На рис. 8 дано сравнение распределения локального коэффициента трения су по поверхности пластины конечной длины с блазиусовским распределением су на отрезке длины Ь

—1/2 1/2

полубесконечной пластины [19] с у ЯеЬ' = 0.332х' , где х — безразмерная продольная

координата,

отсчитываемая от передней кромки. Эта зависимость показана штриховой линией. Сплошные линии соответствуют озееновскому обтеканию, точечные кривые — численному решению уравнений Навье — Стокса [11]. Кривые 1 получены при Re^ = 1, 2 — при Re^ = 10. Видно хорошее согласие между соответствующими кривыми, основное различие наблюдается вблизи концов пластины. С увеличением Re^ распределение Cf приближается к блазиусовскому.

Автор благодарит Ю. Е. Маркачёва и А. Г. Наливайко за техническую помощь. Особая благодарность D. H. Bailey (NASA Ames) за создание и размещение в Интернет пакета программ для вычислений с повышенной точностью и программы-транслятора к нему.

Работа выполнена при финансовой поддержке РФФИ, грант 05-01-00646.

ЛИТЕРАТУРА

1. Tomotika S, Aoi T. The steady flow of a viscous fluid past a sphere and circular cylinder at small Reynolds numbers // Quart. J. Mech. Appl. Math. — 1950. Vol. 3, N 2.

2. Tomotika S, Aoi T. The steady flow of a viscous fluid past an elliptic cylinder and a flat plate at small Reynolds numbers // Quart. J. Mech. Appl. Math. — 1953. Vol. 6, N 3.

3. Hasimoto H. On the flow of a viscous fluid past an inclined elliptic cylinder at small Reynolds numbers // J. Phys. Soc. Jap. — 1953. Vol. 8, N 5.

4. Yamada H. On the slow motion of viscous liquid past a circular cylinder // Rep. Res. Inst. Appl. Mech. Kyushu Univ. — 1954. Vol. 3, N 9.

5. B ourot J.-M. Sur le calcul numerique des champs hydrodynamiques d’Oseen, autour d’une sphere, pour differents nombres de Reynolds // J. Mec. — 1970. Vol. 9. N 4.

6. A bane S. M. Calculation of Oseen flows past a circular cylinder at low Reynolds numbers // Appl. Sci. Res. — 1978. Vol. 34, N 4.

7. Кашеваров А. В. Обтекание сферы и кругового цилиндра в приближении Озеена // Ученые записки ЦАГИ. — 2000. Т. XXXI, № 1 — 2.

8. Dennis S. C. R., Chang G.-Z. Numerical integration of the Navier-Stokes equations for steady two-dimensional flow // Phys. Fluids. — 1969. Vol. 12, N 12, Suppl. II.

9. Меллер Н. А. Обтекание эллиптического цилиндра потоком вязкой жидкости // ЖВММФ. — 1978. Т. 18, № 2.

10. Hudson J. D., Dennis S. C. R. The flow of a viscous incompressible fluid past a normal flat plate at low and intermediate Reynolds numbers: the wake // J. Fluid Mech. — 1985. Vol. 160.

11. Dennis S. C. R., Dunwoody J. The steady flow of a viscous fluid past a flat plate // J. Fluid Mech. — 1966. Vol. 24, N 3.

12. Seebass R., Tamada K., Miyagi T. Oseen flow past a finite flat plate // Phys. Fluids. — 1966. Vol. 9, N 9.

13. Miyagi T. Numerical study of Oseen flow past a perpendicular flat plate // J. Phys. Soc. Jap. — 1968. Vol. 24, N 1.

14. Miyagi T. Standing vortex pair behind a flat plate normal to uniform flow of a viscous fluid // J. Phys. Soc. Jap. — 1978. Vol. 45, N 5.

15. Мак-Лахлан Н. В. Теория и приложения функций Матье. — М.: Изд. иностр. лит. - 1953.

16. Bailey D. H. A multiple precision floating point computation package // Trans. Math. Software. — 1993. Vol. 19, N 3.

17. Таблицы для вычисления функций Матье / Библ. матем. табл., вып. 42. — М.: ВЦ АН СССР. — 1967.

18. Keller H. B., Takami H. Steady two-dimensional viscous flow of an incompressible fluid past a circular cylinder // Phys. Fluids. — 1969. Vol. 12, N 12, Suppl. II.

19. Шлихтинг Г. Теория пограничного слоя. — М.: Наука. — 1974.

Рукопись поступила 1/VII2005 г.