CC BY
58
Вычислительные технологии
Том 15, № 5, 2010
Моделирование начальной стадии отрывного обтекания разомкнутого контура методом дискретных вихрей
Д. Н. Горелов Омский филиал Института математики СО РАН, Россия e-mail: gorelov@ofim.oscsbras.ru
А. И. Говорова Омский государственный университет, Россия e-mail: govorovaai@mail.ru
В нелинейной теории крыла в плоском нестационарном потоке одной из нерешенных проблем является решение задачи в начальной стадии движения крыла. В настоящей работе предложено решение этой проблемы для отрывного обтекания разомкнутого контура в рамках модели идеальной несжимаемой жидкости. Начально-краевая задача сводится к решению системы нелинейных интегродиф-ференциальных соотношений при заданных начальных условиях. Новыми элементами в постановке задачи являются нелинейные дифференциальные соотношения в точках схода вихревых следов с контура и условия непрерывности вихревых слоев, моделирующих контур и вихревые следы, в тех же точках. Введение этих условий позволило получить нелинейное уравнение для определения координат свободных дискретных вихрей, сходящих с кромок контура в рассматриваемые моменты времени. Построена асимптотика решения в окрестности начального момента времени, разработан алгоритм расчета координат свободных дискретных вихрей, скорости схода и интенсивности вихрей, сходящих с контура в вихревые следы. Приведены результаты численного эксперимента для задачи отрывного обтекания пластинки, начавшей движение с постоянной скоростью из состояния покоя.
Ключевые слова: нелинейная начально-краевая задача, отрывное обтекание контура, вихревые следы, метод дискретных вихрей.
Задача нестационарного отрывного обтекания разомкнутого контура решалась многими авторами. Главной ее особенностью является существенная нелинейность течения жидкости вблизи кромок контура в начальные моменты времени, когда начинают формироваться вихревые следы. Для некоторых модельных течений была изучена асимптотика решения в малой окрестности начального момента времени (см., например, [1, 2]). Но основным направлением исследований остается численный эксперимент. В известных алгоритмах решения вихревые следы моделируются системой свободных дискретных вихрей, а задача Коши решается с применением процедуры пошаговой дискретизации по времени. При этом определение координат свободных вихрей, сходящих с кромок контура в начальные моменты времени, остается проблемой с середины XX века. Обычно эти координаты задаются априорно. В настоящей работе предложен алгоритм, который дает решение задачи для начальной стадии отрывного обтекания контура, включая определение координат свободных дискретных вихрей.
1. Постановка задачи
Рассмотрим нелинейную начально-краевую задачу нестационарного отрывного обтекания разомкнутого гладкого контура Ь, начавшего движение с постоянной скоростью в идеальной несжимаемой жидкости, В системе координат 0ху, связанной с контуром Ь, скорость жидкости в бесконечном удалении от контура Ь равна = 0 для £ < 0 и = да я £ > 0. В такой постановке задачи для момен тов времени £ > 0 циркуляция скорости Г(£) вокруг контура Ь меняется с течением времени, в результате чего с кромок контура (точки А, В начинают сходить вихревые следы Ьу1,Ьу2 (рис, 1), Границами области течения являются бесконечно удаленная точка г = оо, г = х + гу, контур Ь и вихревые следы Ьуз- (£), которые эволюционируют с течением времени ¿. Следуя [3], начально-краевая задача для комплексной скорости г)(г,£) может быть сформулирована следующим образом. Вне контура Ь и вихревых следов Ьу= 1, 2, комплексная скорость представима в виде
1
ън
7(в,¿) *
Е
3=1
1
ЪН
7¿а г - гу (а,£)
(1)
0 зВ = 5*2 = 1 интенсивности
где ((з,£) Е Ь; гу(а, ¿) € Ьуз-(£); в, а — дуговые координаты (зА = гШ1(0,£) = г а гад2(0,£) = ¿в); 1 — длина контура Ь; 7 (з,£), 7у (а,£) вихревых слоев, моделирующих контуры Ь,ЬУ3-.
В начальный момент времени £ = 0 комплексная скорость г)(г, 0) = 0, циркуляция скорости Г(0) = 0, а вихревые следы отсутствуют. Граничными условиями являются: затухание возмущенных скоростей в бесконечно удаленной точке (выражение (1) удо-
Ь
г Е Ь)
Тш < е
1
2^
7(5,¿)
*-СМ
21 У —
^ 2пг 3 = 1
¿а
г - гу (а,£)
(2)
и условия непрерывности давления и нормальной составляющей скорости жидкости на Ьуз-, Граничные уеловия па Ьу выполняются при свободном перемещении точки
гуз Е Ьу вместе с жидкостью. Это имеет место, если комплексная координата гу (а, ¿) вихря, сошедшего с точки г*3- контура Ь (г+1 = г^, г*2 = гв) в некоторый момент времени т (0 < т < £), является решением задачи Коши
Ь
у2
Ь
й
— = 2^(0, т) = те [0,г), у(г,0) = 0. (3)
Здесь v(zwj(а,г),г) определяется формулой (1), в которой сингулярные интегралы по контурам Lwj (г) понимаются в смысле главного значения по Коши, К соотношениям (1)-(3) следует добавить условие постоянства суммарной циркуляции скорости по контурам L, Lwl, Lw2 (теорема Кельвина о постоянстве циркуляции скорости по замкнутому жидкому контуру):
£ Ш
Г(г) +(г)
0, Г(г) = у 7(з,г)йз, ^(г)= ] Ywj(а,г)йа. (4)
Ь Ьп, А
Кроме того, будем требовать ограниченность комплексной скорости в кромках контура L (постулат Жуковского) и в концевых точках вихревых следов. Постулат Жуковского выполняется при условиях схода вихревых следов с кромок по касательной к L и непрерывности вихревого слоя в точках схода:
7(з^,г) = Ywj(0,г), 3 = 1, 2, 8*1 = 0, 5*2 = I. (5)
Начально-краевая задача (1)-(5) в общем случае не имеет аналитических решений, В настоящее время задачи такого типа обычно решают путем применения процедуры дискретизации по времени и моделирования вихревых следов системой свободных дискретных вихрей, перемещающихся вместе с жидкостью [4-9], Это позволяет свести исходную нелинейную задачу с неизвестными границами, зависящими от времени, к последовательному решению краевых задач в фиксированных областях для ряда значений времени г1,...,гп.
Переход от одного значения времени к другому требует достаточной гладкости решения в зависимости от г. Это, как правило, имеет место для гп,п> 2. Однако на первых шагах, особенно для г1, исходная начально-краевая задача может быть настолько нелинейной, что ее решение становится проблематичным.
Цель настоящей работы — построение алгоритма решения исходной нелинейной задачи для начальной стадии обтекания контура, включая определение координат свободных дискретных вихрей и основных характеристик отрывного нестационарного течения в точках схода вихревых следов с контура.
2. Общий алгоритм решения
В соответствии с методом дискретных вихрей [4, 7] контур L моделируем системой дискретных вихрей Гт(г), т = 1,..., N расположенных в точках zm Е L, а вихревые следы Lwj (гп), сошедшие с кон тура ^ ^^ ^емя ¿п, моделируем соответственно системой свободных дискретных вихрей Г3ги)к, к = 1,..., п, 3 = 1, 2, которые в момент времени Ьп находятся в точках z^ji>k(tn), Величины Г3шк определяют суммарную интенсивность элемента ДLjUlk (г) вихревого сл еда Lwj•, сошедшего с кон тура ^ ^^ ^емя Дгк = гк — гк-1. Этот элемент деформируется с течением времени, но его суммарная интенсивность остается постоянной. При таком моделировании контура L и вихревых следов формула (1) для комплексной скорости принимает вид
т= 1 ]=1 к=1 wk\ )
Потребуем выполнения условия (2) непротекания жидкости через контур Ь в N +1 контрольных точках 20р € Ь, р =1,..., N + 1, полагая 201 = 2*ь 20^+1 = 2*2, а дискретные вихри Гт разместим между контрольными точками. Тогда для момента времени ¿п получим систему из N + 1 уравнений (р = 1,..., N + 1):
1 N Г (7 ) 1 2 п Г /т{е* (*<»,*„) [г) + у ^^ + у У -Ь*-]} = о. (7)
Условие (4) постоянства суммарной циркуляции скорости по контурам Ь, Ьш1, Ьад2 можно записать в виде
N 2
- Гта(<п-1)] = 0. (8)
т=1 .7=1
Искомыми величинами в системе N + 2 уравнений (7), (8) являются: интенсивности дискретных вихрей Гт(£п), т = 1,..., N на контуре Ь, интенсивности свободных дискретных вихрей Г^п, ] = 1, 2, и координаты всех свободных дискретных вихрей в момент времени ¿п. Но уравнения (7), (8) позволяют определять только интенсивности дискретных вихрей Гт(£п), Г^п при условии, что координаты свободных вихрей Г^п известны,
В настоящее время система уравнений (6)-(8) является основной при решении задач нестационарного обтекания разомкнутого контура методом дискретных вихрей, В известных алгоритмах решения рассматриваемой задачи применяются следующие процедуры, Для момента времени ¿1 задаются координаты 1(^1) свободных вихрей Г^1, после чего решением системы уравнений (7), (8) определяются интенсивности дискретных вихрей Гт(^) на контуре Ь и свободных дискретных вихрей Г^1. Затем по формуле (6) проводится расчет скорости жидкости в точках 1(^1), на основе которого с помощью уравнения (3) вычисляются новые координаты вихрей Г^ в момент времени ¿2. В дальнейшем для ¿п, п > 2, расчет ведется по тому же алгоритму, включая задание координат 2^п(^п) свободных вихрей Г^п, сходящих с кон тура Ь, При этом для каждого момента времени ¿п, п > 2, интенсивности свободных вихрей к = 1,...,п — 1, известны из решения задачи в предыдущие моменты времени, а их координаты в момент времени ¿п определяются решением уравнения (3),
Таким образом, в известных алгоритмах численного решения рассматриваемого класса задач координаты свободных вихрей Г^п в момент времени ¿п не определяются, а задаются априорно.
3. Решение задачи в начальные моменты времени
Для определения координат свободных дискретных вихрей, сходящих с контура в моменты времени ¿п, к формулам (6)-(8) необходимо добавить новые независимые уравнения, в качестве которых выберем соотношения
- 7^(0, = 0, 3 = 1, 2, (9)
Ь
Здесь ад, (¿) — скорость схода вихревого следа Ь^ с Ь,
Условия (5), (9) записаны для непрерывного вихревого слоя, моделирующего контур Ь и вихревые следы Система же уравнений (6)-(8) построена для определения ин-тенеивноетей дискретных вихрей Гт(Ь) и Г,га, определяющих суммарные интенсивности вихревых слоев па соответствующих элементах контуров Ь и Ь,, Поэтому для получения единой системы уравнений для определения интенеивноетей дискретных вихрей и координат г3шк(¿), которые входят в формулу (6), необходимо функции 7(^, ¿), 7,- (0, ¿), Wj (¿) и Гwj (¿) представить через интенсивности дискретных вихрей Гт, Г,га.
Предварительно исследуем асимптотику рассматриваемых функций в малых окрестностях точек (концов контура Ь) и начального момента времени, полагая Ь € [0, ¿1]. Комплексная скорость, индуцируемая элементом вихревого следа ДЬ,1(Ь), определяется интегралом Коши (1), В классе функций, ограниченных на концах разомкнутого контура, плотность этого интеграла (функция 7,-(а, ¿)) должна быть равна нулю на конце контура ДЬ,1(Ь) с порядком убывания [1 — а//,-(¿)]1/2, а € [0,/,-(¿)], где lwj(¿) — длина элемента ДЬ,1(Ь). Это следует, в частности, из формул обращения интеграла типа интеграла Коши [10].
В соответствии с поведением интенсивности вихревого слоя в конце следа функция 7,- (а, ¿) на элементе ДЬ,1(Ь) в первом приближении меняется по параболическому закону:
7wj(а,Ь) = 7wj(0,Ь)[1 — а/^(¿)]1/2, а € [0,lwj(¿)], Ь € [0,^]. (10)
В этом случае циркуляция скорости по элементу ДЬ,1(Ь)
I ^¿(а,^ = Ь^ = аь3ш1. (п)
3
Дифференцируя (11) по t, имеем
drwj 2
dt ~ 3
d d
Lj(t)--fwj(0,t) + jwj(0,t)—lwj(t)
(12)
Предположим, что для Ь € [0,Ь1] длина lwj• (¿) элемента ДЬ,1(Ь) изменяется только за
Ь
dl
= Wj, Wj = const. (13)
Из (9)—(13) следует, что для t £ [0, 11]
lW3(t) = w3t, ^ = 2 (14)
Сравнивая (14) с (9), придем к уравнению 2d^wj/jwj = dt/t, решением которого является 7^(0, t) = СлД, С = const.
Полученные соотношения позволяют сделать вывод, что при параболическом законе (9) изменения интенсивности вихревого слоя на элементе ALW1(t^ и Wj (t) = const для t £ [0, ti]
lwj = Wjt, 7^(0, t) = С Vi, rw3 = - Cw3t3/2, = 2 (15)
Выражения (15) определяют асимптотику соответствующих функций в окрестности начального момента времени Ь € [0,^],
Тем же путем можно исследовать асимптотику решения вблизи кромок контура Ь для моментов времени Ь € [¿п-1,£п], п > 2. В этом случае интенсивность вихревого следа (а, ¿) на элементе АЬ2П(Ь) мало отличается от значения 7^3 (0, что позволяет полагать
7^3(а,Ь) = 7^3(0,Ь), а € [0,/у(Ь)], Ь € [¿п-1,^п]. Тогда асимптотические выражения (15) переходят в следующее:
и = еош^ /23 = и(Ь - Ьп-1), 72(0,Ь)= 7(з*3,Ь), Г^ = 7^3(0,Ь)и(Ь - ¿^-1),
^Г •
= ¿е[*га_ь*га]. (16)
Вернемся к условиям (5), (9), раеемотривая их в момент времени Из (15) следует, что /^ = 3/2 Г^/АЬ, Скорость схода и (Ь1) вихрей в след определяется формулой (6) как проекция вектора относительной скорости жидкости на касательные к контуру Ь в точках ¿а, ¿в- Условие непрерывности вихревого слоя (5) в кромках контура Ь позволяет заменить 7^3(0,Ь1) на функцию 7(3*3,Ь1), которую можно выразить через дискретные вихри Гт(Ь1), применяя подходящую интерполяционную формулу, В частности, при достаточной гладкости функции 7(3) вблизи концов контура Ь можно ограничиться формулами 7^1(0,^1) = 7(3*1,^1) = Г^^/А, 7^2(0,^1) = 7(3*2,^1) = Гм(^)/А, где А = //Ж, В результате все величины в (9) оказываются зависящими только от интенсивности дискретных вихрей Г1(Ь1),..., Г^ и Г^1, Г^1. При заданном положении свободных вихрей Г^1, Г^1 интенсивности всех этих вихрей определяются решением уравнений (7), (8), но при этом условия (13) будут выполняться только при некоторых значениях координат являющихся искомыми величинами в рассматриваемой
задаче.
Построим алгоритм расчета координат свободного дискретного вихря Г^, при которых могут быть выполнены все условия задачи (5)-(9). Представим искомые координаты в виде
31 = ¿*3 + (—е Сад?,
= + ^ = (¿13 + ¿¿2,-)А, э = 1, 2, А = //Ж, (17)
где — значения угла в в концах контура Ь, а безразмерные величины ¿3, ¿23- определяют комплексные координаты вихря Г^ в системе коордипат , п23 (рис, 2), Из (17) следует, что каждой паре величин ¿у, ¿3 соответствуют свои комплексные координаты 1 свободных вихрей, сходящих с контура Ь за время АЬ, Подставляя эти
Рис. 2. Моделирование элемента ДЬ^^х) вихревого следа дискретным вихрем
координаты в уравнения (7) и решая данные уравнения совместно с (8), определим интенсивности дискретных вихрей ri(ti),..., rN(ti) и Г^, после чего с помощью (6) вычислим скорости схода Wj (ti) с контура L, Полученные результаты позволяют вычислить функции
fj(Sij, ¿2j) = ^rWi - Wj(ti)7(s*j,ti)Aii, j = 1, 2, ß = 3/2. (18)
При произвольных значениях öij, ö2j функции fj (öij ,ö2j) = 0 и определяют в принятом приближении невязку соотношений (9), Пусть пара величин öij, ö2j обращает функцию fj в ноль. Тогда при соответствующих им координатах z3wi свободных дискретных вихрей задача отрывного обтекания контура решается с выполнением всех поставленных условий (5)-(9).
Таким образом, задача определения координат свободных вихрей, сходящих с контура L за время At = ii, свелась к расчету корней нелинейного уравнения (18). Аналогичным образом можно построить решение задачи определения комплексных координат zWn, zWn для любого момента времени tn, n > 2, но при этом в соответствии с (16) коэффициент ß в формуле (18) должен быть равен ß = 1.
4. Численный эксперимент
Предложенный алгоритм был применен для решения задачи отрывного нестационарного обтекания пластинки, начавшей движение с постоянной скоростью из состояния покоя. Пластинка была установлена под углом в = 90° к вектору скор ости v^. Шаг по времени выбирался равным At = tn — tn-i = T/N, T = //|v^|, N = 20.
Нелинейное уравнение fj (öij ,ö2j) = 0 решалось численно путем перебора переменных öij, ö2j в некоторой заданной области с учетом градиента функции fj по этим переменным. Расчет показал, что уравнение fj (öij ,ö2j) = 0 имеет бесконечное множество корней öij, ö2j, определяющих координаты свободных вихрей, сходящих с контура L за At
таты расчета безразмерных координат öij ,ö2j свободных вихрей, сходящих с контура в моменты времени ti, t2, t3, и соответствующие им скорости схода вихревых следов с
ti
ные значения ö2j оказываются положительными, тогда как в последующие моменты времени координата ö2j может иметь отрицательные значения. На рисунке представлены данные только для ö2j > 0 так как зпачения ö2j < 0 в рассматриваемой задаче отрывного обтекания пластинки не имеют физического смысла.
Таким образом, нелинейная задача (6)-(9) методом дискретных вихрей имеет бесконечное множество решений. В связи с этим возникает проблема выделения из бесконечного множества корней öij ,ö2j уравнения fj (öij ,ö2j) = 0 той пары значений öij ,ö2j, для которой решение рассматриваемой задачи имеет физический смысл. В общем случае каждому значению öij при öij G (0, соответствуют два значения ö2j, то при öij = оба эти значения становятся равными ö2j = ¿2*. Иначе говоря, при
öij = ¿i*, ¿2j = ¿2* (19)
решение задачи (6)-(9) становится единственным и имеющим физический смысл. На рис. 3 решения, соответствующие координатам (19), показаны точками.
а б
Рис. 3. Результаты расчета безразмерных координат ¿13 ,¿23 свободных дискретных вихрей, сходящих с контура в моменты времени ¿1, ¿2, ¿3 (я)5 и соответствующие им скорости схода следов с контура (б)
На рис. 4 представлены результаты расчета координат п свободных дискретных вихрей (точки), сходящих с верхней кромки пластинки (э = 2), для моментов времени Ь1,..., Ь5, Приведенные данные позволяют сделать вывод, что в первые моменты времени, когда след только зарождается, основное влияние на форму следа оказывает набегающий поток. Но с течением времени след растет и начинает жить самостоятельно, что проявляется в его интенсивной эволюции и в тенденции к скручиванию. Расчет показал, что начальный участок вихревого следа, примыкающий к кромке контура, практически не деформируется. Это обстоятельство позволяет упростить алгоритм расчета неетаци-
Рис. 4. Развитие вихревого следа, сходящего с верхней кромки пластинки, в начальные моменты времени ¿1, ¿2, ¿3, ¿4, ¿5
онарного отрывного обтекания контура, полагая координаты свободных вихрей Г2п, сходящих с контура в момент времени Ьп, п > 2, заданными (из раечета для Ь = Ь1).
Как уже отмечалось, в известных алгоритмах координаты ¿¿п свободных вихрей Г2п, ^ящих с контура в каждый момент времени Ьп, задаются априорно. Наиболее часто сходящие вихри размещают по касательной к контуру на расстоянии А/2, А = от кромки. В представленном примере этому положению на рис. 4 соответствует кружочек на оси ординат. При таком задании координат сходящих вихрей соотношения (9) не выполняются, а эволюция вихревых следов будет отличаться от данных на рис. 4. Оценить точность решения рассматриваемой начально-краевой задачи при таком задании координат ¿¿п трудно. Для этого нужно сравнивать не только кинематическую картину течения, но и гидродинамические реакции на контур. В настоящей работе описана только кинематика отрывного нестационарного течения. Расчет гидродинамических реакций на контур авторы планируют провести в следующей работе.
Заключение
Разработан алгоритм решения начальной стадии отрывного обтекания разомкнутого контура методом дискретных вихрей в рамках нелинейной теории крыла в нестационарном потоке, позволяющий решать задачу в полном объеме, включая расчет координат свободных дискретных вихрей, сходящих с контура. Традиционная постановка задачи была дополнена требованием выполнения в кромках контура дифференциальных соотношений, связывающих цикуляцию скорости вокруг вихревого следа, интенсивность вихрей и скорость их схода с контура. С помощью этих соотношений построена асимптотика решения в окрестности начального момента времени, определены координаты свободных дискретных вихрей и локальные характеристики течения. Установлено, что рассмотренная нелинейная задача допускает бесконечное множество решений. Предложен способ выделения из этого множества единственного решения.
Список литературы
[1] Никольский A.A. О второй форме движения идеальной жидкости около обтекаемого тела (исследование отрывных вихревых потоков) // Докл. АН СССР. 1957. Т. 116, № 2. С. 193-196.
[2] Бетяев С.К. Эволюция вихревых пелен //Динамика сплошной среды со свободными поверхностями: Сб. науч. тр. Чебоксары: Чувашский гос. ун-т, 1980. С. 27-38.
[3] Горелов Д.Н. К постановке нелинейной начально-краевой задачи нестационарного отрывного обтекания профиля // ПМТФ. 2007. Т. 48, № 2. С. 48-56.
[4] Сарпкайя Т. Вычислительные методы вихрей. Фримановская лекция (1988) // Тр. Американского общества инж.-мех. Машиностроение. Сер. А. 1989. № 10. С. 1-60.
[5] Giesing J.P. Nonlinear two-dimensional unsteady potentional flow with lift //J. Aircraft. 1968. Vol. 5, No. 2. P. 135-143.
[6] Ильичев К.П., Постоловский С.Н. Расчетное исследование нестационарного отрывного обтекания тел плоским потоком невязкой жидкости // Изв. АН СССР. МЖГ. 1972. № 2. С. 72-82.
[7] Белоцерковский С.М., Ништ M.II. Отрывное и безотрывное обтекание тонких крыльев идеальной жидкостью. М.: Наука, 1978. 352 с.
[8] Головкин В.А. Нелинейная задача о неустановившемся обтекании произвольного профиля со свободно деформирующимся вихревым следом // Уч. зап. ЦАГИ. 1972. Т. 3, № 3. С. 1-11.
[9] Молчанов В.Ф. Некоторые вопросы расчета течений с тангенциальными разрывами// Уч. зап. ЦАГИ. 1975. Т. 6, № 4. С. 1-11.
[10] Гахов Ф.Д. Краевые задачи. М.: Физматгиз, 1963. 640 с.
Поступила в редакцию 10 июня 2009 г., с доработки 18 января 2010 г.
