CC BY
38
Наука к Образование
МГТУ им. Н.Э. Баумана
Сетевое научное издание
ISSN 1994-0408
Наука и Образование. МГТУ им. Н.Э. Баумана. Электрон. журн. 2015. №2. С. 234-249.
DOI: 10.7463/0215.0756954
Представлена в редакцию: 19.02.2015 © МГТУ им. Н.Э. Баумана
УДК 532.5+519.6
О точности расчетных схем вихревых методов при моделировании обтекания профилей с угловой точкой
Кузьмина К. С.1, Марчевский И. К.1'*, Морева В. С.1
1 Россия, МГТУ им. Н. Э. Баумана
Рассмотрена задача о моделировании двумерного обтекания профиля несжимаемой средой с использованием метода вихревых элементов. Проведено исследование точности расчетных схем, основанных на разных подходах к удовлетворению граничного условия, применительно к моделированию жесткого неподвижного профиля с угловой точкой. Отличие данных подходов состоит в способе удовлетворения граничного условия, которое обеспечивается равенством нулю либо нормальных, либо касательных компонент скорости на профиле. Исследовано влияние равномерности разбиения профиля на панели на точность расчетных схем. В случае, когда разбиение профиля является близким к равномерному, расчетные схемы, основанные на обоих подходах, позволяют получить качественно правильное решение с близкой достаточно высокой точностью. В случае существенного различия длин панелей, примыкающих к угловой точке, метод с нормальными компонентами, в отличие от метода с касательными компонентами скорости, не позволяет получить качественно правильное решение.
Ключевые слова: сходимость; интегральное уравнение; метод вихревых элементов; расчетная схема; особенность решения
При решении двумерных задач аэрогидродинамики и аэрогидроупругости в случаях, когда скорость относительного движения обтекаемого профиля и среды является малой и влиянием сжимаемости среды можно пренебречь, весьма удобными для практического использования могут быть бессеточные лагранжевы вихревые методы [1, 2, 3, 4, 5]. Они особенно эффективны при моделировании внешнего обтекания профилей, когда область течения является неограниченной. Также существенным преимуществом вихревых методов является то, что расчет обтекания подвижных и/или деформируемых профилей не является принципиально более сложным по сравнению со случаем неподвижного жесткого профиля.
Введение
Течение несжимаемой жидкости в неограниченной области описывается уравнениями
дУ - V»
У-У = 0, д- + (V -У)У = vДV -Ур, (1)
дЬ р
где V — скорость; р'— давление; р — плотность; V — кинематическая вязкость среды.
Начальное распределение скоростей среды предполагается известным; на бесконечности задается условие затухания возмущений, а на поверхности обтекаемого профиля — условие прилипания (либо непротекания в случае V = 0):
у ^ У \?€К = (или У \?€К ■ п = -и). (2)
Здесь УК — скорость движения точек поверхности профиля, которая считается известной; Уп — орт нормали к поверхности профиля.
Решение задач вихревыми методами предполагает моделирование профиля распределенным по его поверхности вихревым слоем и слоем источников. Известны также и другие подходы, предполагающие, в частности, наличие завихренности внутри профиля [2]. Слой источников является присоединенным, т. е. связанным с профилем; его рассмотрение необходимо в случае подвижного (деформируемого) профиля. Вихревой слой в общем случае состоит из двух частей — присоединенного и свободного слоев; свободная завихренность генерируется на поверхности профиля и сходит в поток, образуя вихревой след вблизи и позади профиля. Интенсивности вихревых слоев и слоя источников определяются из граничного условия прилипания либо непротекания на поверхности профиля (соответственно, при моделировании течения вязкой либо идеальной среды). Известно несколько подходов к построению расчетных схем вихревых методов, отличающихся способом удовлетворения граничного условия [6, 7, 8].
Таким образом, расчет обтекания профиля включает в себя решение двух основных задач — определение интенсивностей вихревого слоя и слоя источников, а также моделирование эволюции вихревого следа в области течения. Корректное вычисление интенсивностей слоев имеет определяющее значение для правильности решения задачи в целом. Целью настоящей работы является исследование точности соответствующих расчетных схем в случае моделирования обтекания профилей с угловой точкой.
1. Краткое описание расчетных схем метода вихревых элементов
Ограничимся рассмотрением наиболее простого случая, когда обтекаемый профиль является жестким и неподвижным. При этом присоединенный слой источников будет отсутствовать. Границу обтекаемого профиля будем считать кусочно-гладкой замкнутой кривой К, орты внешней нормали и касательной к ней (на гладких участках) обозначим за п и V соответственно.
В качестве модельной задачи рассмотрим ситуацию, когда область течения является неограниченной, завихренность в ней отсутствует, а набегающий поток, имеющий скорость Т/^, возмущается лишь профилем. Такая постановка задачи соответствует мгно-
венному старту течения среды. Задачи расчета интенсивности вихревого слоя, которые необходимо решать на каждом шаге моделирования нестационарного обтекания профиля, отличаются от описанной лишь тем, что в области течения имеется известное распределение завихренности, которое будет оказывать влияние на правую часть соответствующих интегральных уравнений.
«Классический» подход к удовлетворению граничного условия состоит в обеспечении равенства нулю нормальной компоненты скорости среды на поверхности профиля и приводит к необходимости решения интегрального уравнения [1,6]
1 х (г - Г) ■ й(О) 7(0 ^ = -Ко ■ п(г), г е К, (3)
2п / \ 1г - б I2
к у 1 ^ 1
где к — единичный вектор, перпендикулярный плоскости течения, выбираемый таким образом, что к х п(г ) = Г (г ); 7 (б ) — искомая интенсивность вихревого слоя на поверхности профиля.
Данное уравнение является сингулярным с ядром Гильберта, интеграл понимается в смысле главного значения по Коши. Его решение неединственно, для выделения нужного решения используется дополнительное условие
Y(V) ^ = 0. (4)
K
Реализации вихревых методов, основанные на решении уравнений (3) и (4), будем обозначать аббревиатурой НМВЭ — метод вихревых элементов с нормальными компонентами скорости.
Альтернативный подход к разработке расчетных схем вихревых методов основан на равенстве нулю предельного значения касательной компоненты скорости со стороны профиля [6, 7, 8]
V- (r) ■ V(r ) = V(r ) ■ V(r) - = 0, r e K, (5)
которое приводит к необходимости решения интегрального уравнения вида
i /(к -■ V(V))y(V) dl« - ^ = -V- ■ V(V), r e K. (6)
Как и в предыдущем случае, данное уравнение имеет неединственное решение. Для выделения единственного решения используется то же самое условие (4). Расчетные схемы вихревых методов, предполагающие решение уравнений (6) и (4), обозначим аббревиатурой КМВЭ — метод вихревых элементов с касательными компонентами скорости. Отметим, что уравнения (3) и (6) эквивалентны [6].
Ядро интегрального уравнения (6) является ограниченным, если K является гладкой кривой класса C2, а также при моделировании обтекания профилей с острой кромкой, являющейся точкой возврата (например, семейство профилей Жуковского). В этих случаях решение уравнения (6) будет ограниченной функцией.
Расчетные схемы, предполагающие численное решение уравнения (6), показали существенно более высокую точность по сравнению со схемами, основанными на уравнении (3), при решении различных модельных задач [7, 8, 9, 10]. В задачах, рассмотренных в [7, 8, 9] осуществлялся расчет интенсивности вихревого слоя на эллиптических профилях и профилях Жуковского как в описанной выше постановке, так и при наличии завихренности в области течения. Выбор этих профилей был связан с двумя факторами: как указано выше, для таких профилей ядро уравнения (6) является ограниченным, и для них несложно построить точное аналитическое решение задачи, используя аппарат теории функций комплексного переменного [11, 12].
В то же время известно [1], что при наличии на профиле угловой точки с внешним углом а > п, интенсивность вихревого слоя в окрестности этой точки может иметь слабую (интегрируемую) особенность вида рп/а-1. Отсюда следует, что интенсивность вихревого слоя возрастает при приближении к острой кромке и вблизи нее в поток могут сходить вихри с достаточно большой циркуляцией, что подтверждается результатами практических расчетов. Точность вычисления этих циркуляций во многом определяет правильность моделирования вихревого следа.
В настоящей работе рассмотрим модельную задачу о расчете обтекания профиля с угловой точкой и исследуем точность определения циркуляций вихревых элементов, генерируемых вблизи угловой точки.
Рассмотрим обтекание симметричного профиля, имеющего длину Ь и характеризуемого углами а и в как показано на рис. 1.
2. Постановка задачи
V.
ь
Рис. 1. Профиль с угловыми точками
Я
вт(в/2) — сое а'
задняя часть образована отрезками прямых. Набегающий поток является горизонтальным. При этом искомая интенсивность вихревого слоя будет обращаться в бесконечность в верхней и нижней точках профиля (предполагается, что а > п) и при приближении к ним, как следует из [1], 7 ~ Срж/а-1, где р — расстояние до угловой точки, С — некоторая константа.
Построение расчетной схемы в вихревых методах предполагает замену исходного профиля ломаной, звенья которой называют панелями, при этом на каждой панели генерируется один вихревой элемент. Ясно, что точность решения зависит от количества панелей; при этом равномерность разбиения может существенно влиять на результат. В частности, в методе НМВЭ длины соседних панелей не должны существенно различаться. Это связано с особенностями используемой квадратурной формулы, позволяющей выделять главное значение в смысле Коши сингулярного интеграла в уравнении (3). В методах КМВЭ на гладких участках профиля разбиение на панели может производиться достаточно произвольно, поскольку подынтегральная функция в уравнении (6) является ограниченной. При этом использование схем КМВЭ для эллиптических профилей позволяет повысить точность определения циркуляций вихревых элементов примерно на порядок (при количестве панелей на профиле N « 200), а для профилей Жуковского — до трех порядков в зависимости от параметров профиля и особенностей используемых расчетных схем.
Рассмотрим различные способы аппроксимации профиля панелями и исследуем для них точность расчетных схем НМВЭ и КМВЭ. При этом будем анализировать решение (интенсивность вихревого слоя) на прямолинейном участке профиля при приближении к верхней угловой точке. Для простоты примем, что дуговая часть профиля разбивается на панели одинаковой длины, и рассмотрим различные способы аппроксимации панелями прямолинейного участка:
• панели имеют равные длины, близкие к длинам панелей дуговой части;
• панели имеют равные длины, существенно отличающиеся от длин панелей дуговой части;
• панели имеют разные длины, образующие, к примеру, геометрическую прогрессию, и уменьшаются при приближении к верхней угловой точке.
3. Результаты расчетов для профиля с равномерным разбиением на панели
Подробное описание расчетных схем, относящихся к методам НМВЭ и КМВЭ, приведено в работах [7, 8, 13, 14], там же предложена методика оценки их точности. Наиболее очевидный путь сравнения вычисленного значения интенсивности вихревого слоя с точным решением является некорректным, особенно применительно к профилю с угловой точкой, при приближении к которой данная величина должна неограниченно возрастать. Физический смысл имеет циркуляция вихревого элемента, который сходит в поток и впоследствии участвует в формировании вихревого следа за профилем. Величина циркуляции равна сум-
марной интенсивности вихревого слоя на некотором участке профиля (на панели в схемах КМВЭ и на половинах двух смежных панелей в НМВЭ).
В данном случае точное аналитическое решение задачи неизвестно, известна лишь его асимптотика при приближении к угловой точке, поэтому воспользуемся следующим приемом: для конкретного профиля с фиксированными параметрами Ь, а, в и заданной скоростью набегающего потока У^ численно решим поставленную задачу для различного количества панелей N. При этом будем использовать обе схемы НМВЭ и КМВЭ, разбиение профиля на панели будем осуществлять равномерно.
Считая, что точное решение вблизи угловой точки имеет вид 7*(р) = С рп/а-1, вычислим суммарную завихренность на пр ближайших к угловой точке панелях:
Г* = / 7*(р) ¿р = I Срп/а-1 ¿р = ^(рЙ - рП/а), г =1, 2, ..., пр. (7)
Здесь рг — расстояние от угловой точки до начала г-й панели.
Метод КМВЭ позволяет определять величины — средние по панелям интенсивности вихревого слоя, поэтому на основе найденных выше значений Г* вычислим величину
7* = Г* = ас рй - р1/а (8)
' рг+1 - р% п р%+1 - р%
Далее образуем функцию Ф(С) по следующему правилу
Пр
Ф(С) = Е(7г - 7*)2, (9)
г=1
где суммирование ведется по пр ближайшим к угловой точке панелям; в настоящей работе принято пр = 5. Значение константы С определяется из условия минимизации функции Ф(С).
Для метода НМВЭ методика в целом сохраняется с той лишь разницей, что в расчете вычисляется величина 7», соответствующая среднему значению интенсивности вихревого слоя на половинах смежных панелей, имеющих номера г и г + 1. Соответственно, интегрирование асимптотического решения нужно проводить по отрезкам с концами в точках
Рг + Рг+1 Рг+1 + Рг+2
-и-.
2 2
Результаты определения величины коэффициента С при различном значении числа панелей на дуговой части профиля и соответствующем значении числа панелей на всем профиле представлены в табл. 1.
Видно, что с увеличением числа панелей на профиле значение коэффициента С «стабилизируется», при этом наблюдаемое различие в оценках для методов НМВЭ и КМВЭ не превышает 1 %. Эту разницу можно объяснить описанными выше особенностями методов, а также тем, что точное решение представляет собой сумму
7ехас1(р) = Ср ^-1 + %), (10)
где к(р) — ограниченная функция, наличие которой никак не учитывалось при построении оценки для С.
Вычисленные значения коэффициента С
Число панелей на половине дуговой части Общее число панелей Снмвэ Скмвэ
20 180 0,8081 0,8119
40 360 0,8016 0,8072
60 538 0,7991 0,8058
80 718 0,7983 0,8052
100 896 0,7977 0,8049
150 1344 0,7972 0,8046
200 1792 0,7970 0,8045
300 2688 0,7970 0,8045
Отметим, что при моделировании обтекания профиля в форме ромба, задняя часть которого совпадает с задней частью профиля, изображенного на рис. 1, а передняя является ее симметричным отражением, особенность решения в верхней угловой точке имеет такой же характер, а значения коэффициента С при ней, определяемые методами НМВЭ и КМВЭ по описанной выше методике, различаются менее, чем на 0.1 %.
На рис. 2 для случая N = 896 приведены графики, на которых точками изображены вычисленные значения средних величин интенсивности вихревого слоя на нескольких ближайших к угловой точке участках профиля, а сплошной линией соединены значения величины 7*, найденные на основе полученной оценки коэффициента С. Визуальное различие графиков объясняется тем, что участки профиля, для которых производится расчет в методах НМВЭ и КМВЭ, разные, о чем было упомянуто выше.
Yi
4.5 4.0 3.5 3.0 2.5
2
4
6
8
10
2
4
6
8
10
а б
Рис. 2. Вычисленные значения (точки) средней интенсивности вихревого слоя вблизи угловой точки
в сравнении с асимптотическим решением (линия) для методов НМВЭ (а) и КМВЭ (б)
Y
Средние величины относительных отклонений вычисленных значений средней интенсивности вихревого слоя от асимптотического решения приведены в табл. 2. Осреднение производилось по пр = 5 ближайшим к угловой точке участкам профиля, для которых производился расчет. Напомним, что все расчеты проведены для случая равномерного разбиения профиля на панели.
Отклонения расчетных значений средней интенсивности вихревого слоя
от асимптотического решения
Число панелей 180 360 538 718 896 1344 1792 2688
¿НМВЭ, % ¿КМВЭ, % 2,41 2,12 1,56 1,36 1,26 1,12 1,12 1,00 1,04 0,94 0,95 0,80 0,92 0,85 0,90 0,83
Анализ полученных результатов позволяет сделать вывод о том, что оба метода НМВЭ и КМВЭ позволяют с достаточно высокой точностью рассчитать интенсивность вихревого слоя на профиле с угловой точкой. Данный вывод, однако, справедлив лишь для ситуаций, близких к рассмотренной модельной задаче, когда вблизи обтекаемого профиля в потоке завихренность отсутствует. Исследование точности различных расчетных схем при наличии в потоке вблизи профиля вихревых элементов представляет собой отдельную нетривиальную задачу, частично рассмотренную (для случая кругового профиля) в [14].
4. Результаты расчетов для профиля с неравномерным разбиением на панели
В ряде случаев обеспечивать равенство длин всех панелей на обтекаемом профиле на практике бывает неудобно; в частности, профиль может быть «сконструирован» из отдельных кривых, для которых разбиение на панели уже построено. То же относится к моделированию обтекания деформируемых профилей, когда длины панелей на профиле естественным образом меняются в процессе расчета.
Для описанной выше модельной задачи рассмотрим сначала ситуацию, когда и дуговая, и прямолинейная части профиля (см. рис. 1) разбиты на панели равномерно, но длины этих панелей различны. Примем, что дуговая часть профиля разбита на 200 панелей, а количество панелей на прямолинейных частях будем варьировать, указывая при этом отношение к > 1 длин панелей, примыкающих к угловой точке. В качестве «эталонного» примем асимптотическое решение с ранее определенным значением коэффициента С.
В табл. 3 приведены средние величины относительных отклонений вычисленных значений средней интенсивности вихревого слоя от асимптотического решения при вышеопи-
Т а б л и ц а 3
Отклонения расчетных значений средней интенсивности вихревого слоя от асимптотического решения при различных длинах панелей на дуговой и прямолинейной частях профиля
Отношение длин панелей к 1,10 1,25 1,33 1,50 2,00 3,00 5,00 10,00
Общее число панелей 966 1070 1128 1244 1592 2288 3680 7158
¿НМВЭ, % 1,37 2,74 3,49 4,93 8,95 15,92 27,21 47,65
¿КМВЭ, % 0,98 1,06 1,10 1,19 1,45 1,92 2,53 3,13
санном способе разбиения профиля. Осреднение производилось по пр = 5 ближайшим к угловой точке участкам профиля, для которых производился расчет.
Видно, что погрешность обеих расчетных схем возрастает при увеличении отношения длин панелей, примыкающих к угловой точке, однако при значительной величине этого отношения схема КМВЭ оказывается на порядок более точной. Для наглядности результаты расчета по обеим схемам для значения к = 10.0 показаны на рис. 3.
У;
У;
10 8 6 4 2
• •
4 6 8 10
а
10 8 6 4 2
4 6
б
8 10
Рис. 3. Вычисленные значения (точки) средней интенсивности вихревого слоя вблизи угловой точки при к =10 в сравнении с асимптотическим решением (линия) для методов НМВЭ (а) и КМВЭ (б)
2
2
Расчеты показывают, что равномерность разбиения прямолинейной части профиля на панели не имеет принципиального значения. В частности, если разбиение осуществляется неравномерно, а длины панелей образуют геометрическую прогрессию со знаменателем д = 0.95 ... 0.99, что соответствует сгущению при приближении к верхней угловой точке, то при сохранении постоянного числа панелей на дуговой части и значения к отношения длин панелей, примыкающих к угловой точке, такими же, какие были приняты в табл. 3, результат меняется незначительно. Отклонение от асимптотического значения для схемы КМВЭ практически не изменяется, а для схемы НМВЭ несколько увеличивается, достигая величины в 54 % при д = 0.95, к = 10; в этом случае общее количество панелей на профиле составляет лишь 400 (200 на дуговой части, 200 на двух прямолинейных частях), а не 7158, как при равномерном разбиении всего профиля (см. табл. 3).
При еще более неравномерном разбиении профиля, когда к > 10, относительная погрешность метода НМВЭ продолжает нарастать, и уже можно говорить о получении качественно неверного решения, тогда как погрешность схемы КМВЭ даже в этом случае остается не слишком большой и удается получить качественно правильное решение. Пример числовых значений для разбиения дуговой части профиля на 200 панелей, знаменателя прогрессии д = 0.98, различного количества панелей на прямолинейной части и соответствующих значений отношения к приведен в табл. 4.
Таблица 4
Отклонения расчетных значений средней интенсивности вихревого слоя от асимптотического решения при больших значениях отношения к длин панелей, примыкающих к угловой точке
Число панелей на прямолинейной части 450 500 550 600
Отношение длин панелей к 13,1 21,9 36,3 60,3
Общее число панелей 650 700 750 800
¿НМВЭ, % 69,6 82,0 106,8 135,1
¿КМВЭ, % 3,7 5,6 8,3 11,5
5. Профиль с внешним углом, меньшим развернутого
В заключение рассмотрим профиль, имеющий угловую точку, внешний угол при которой а < п (рис. 4). Алгоритм построения такого профиля остается прежним.
Уоо
Рис. 4. Профиль с угловыми точками
В [1] указано, что в таких угловык точках интенсивность вихревого слоя не имеет особенности. При равномерном разбиении профиля, когда длины всех панелей примерно одинаковы, оба метода НМВЭ и КМВЭ позволяют получить качественно правильное решение. При неравномерном разбиении, соответствующем небольшим значениям параметра к < 2 (в обозначениях предыдущего раздела), решение также остается качественно правильным и практически не изменяется количественно.
Большая неравномерность разбиения профиля приводит к существеной погрешности решения, получаемого методом НМВЭ, тогда как решение, получаемое по схеме КМВЭ, остается близким к исходному. На рис. 5 приведены зависимости средней интенсивности вихревого слоя от расстояния до угловой точки, измеряемого вдоль прямолинейной части профиля. Сплошная линия соответствует решению, полученному при весьма подробном равномерном разбиении профиля на N & 4 000 панелей, точки — неравномерному разбиению; при этом значение N и соответствующее ему значение к варьировалось, расчет проводился по схемам НМВЭ и КМВЭ. Дуговая часть профиля во всех случаях разбивалось равномерно на 200 панелей, длины панелей прямолинейной части образовывали геометрическую прогрессию со знаменателем q = 0.98 (панели измельчались при приближении к угловой точке).
У;
0.7 0.6 0.5 0.4 0.3 0.2 0.1
0.01
0.02
0.03
У;
0.01
0.02
0.03
У;
0.7 0.6 0.5 0.4 0.3 0.2 0.1
• •
0.005 0.010 0.015 0.020
в
У\
0.004 0.008 д
0.012
У\
0.5 0.4 0.3 0.2 0.1
0.005 0.010 0.015 0.020
г
0.004
0.008
0.012
Рис. 5. Вычисленные значения (точки) средней интенсивности вихревого слоя вблизи угловой точки при N = 440, к ^ 2.81 (а, б); N = 500, к ^ 5.39 (в, г); N = 600, к ^ 15.27 (д, е) в сравнении с решением, полученным при равномерном разбиении (линия). Расчеты выполнены по схемам НМВЭ (а, в, д) и КМВЭ (б, г, е)
б
а
У
Р
Р
Р
е
Заключение
Результаты исследования модельной задачи о расчете обтекания профиля с угловой точкой в отсутствие завихренности в области течения позволяют утверждать, что схемы метода вихревых элементов с нормальными и касательными компонентами скорости на профиле (НМВЭ и КМВЭ соответственно) обеспечивают близкую точность при равномерном разбиении профиля панелями. При неравномерном разбиении, когда примыкающие к угловой точке панели существенно различаются по длине, схема КМВЭ обеспечивает значительно
большую точность по сравнению с НМВЭ. В тех случаях, когда отношение длин этих панелей превышает 10, получить качественно правильное решение методом НМВЭ не удается, в то время как при использовании КМВЭ решение остается качественно верным и имеет приемлемую для многих практических целей погрешность.
Сходные результаты получены как для профиля с внешней угловой точкой, вблизи которой точное решение задачи имеет особенность известного характера, так и для профиля с внутренней угловой точки, в которой точное решение особенности не имеет.
Работа выполнена при частичной финансовой поддержке грантов Президента РФ для государственной поддержки молодых российских ученых — кандидатов наук (проекты МК-3705.2014.8 и МК-5357.2015.8).
Список литературы
1. Лифанов И.К. Метод сингулярных интегральных уравнений и численный эксперимент. М.: ТОО «Янус», 1995. 521 с.
2. Андронов П.Р., Гувернюк С.В., Дынникова Г.Я. Вихревые методы расчета нестационарных гидродинамических нагрузок. М.: Изд-во МГУ, 2006. 184 с.
3. Головкин М.А., Головкин В.А., Калявкин В.М. Вопросы вихревой гидромеханики. М.: ФИЗМАТЛИТ, 2009. 264 с.
4. Cottet G.-H., Koumoutsakos P.D. Vortex Methods: Theory and Practice. Cambridge University Press, 2008. 328 p.
5. Lewis R.I. Vortex Element Methods For Fluid Dynamic Analysis Of Engineering Systems. Cambridge University Press, 2005. 592 p.
6. Kempka S.N., Glass M.W., Peery J.S., Strickland J.H. Accuracy Considerations for Implementing Velocity Boundary Conditions in Vorticity Formulations. SANDIA Report SAND96-0583 UC-700. Sandia National Laboratories, Albuquerque, N.M., March 1996. 52 p.
7. Moreva V.S., Marchevsky I.K. Vortex element method for 2D flow simulation with tangent velocity components on airfoil surface // ECCOMAS 2012 6th European Congress on Computational Methods in Applied Sciences and Engineering: Book of proceedings. Vienna, 2012. 14 p.
8. Kuzmina K.S., Marchevsky I.K. On Numerical Schemes in 2D Vortex Element Method for Flow Simulation Around Moving and Deformable Airfoils // Advanced Problems in Mechanics (APM2014): Proc. oftheXLII Summer School-Conference. St.-Petersburg, 2014. P. 335-344.
9. Макарова М.Е. Расчет стационарного безотрывного обтекания профиля потоком идеальной несжимаемой среды // Вестник МГТУ им. Н.Э. Баумана. Сер. Естественные науки. 2011. Спец. вып. Прикладная математика. C. 63-74.
10. Макарова М.Е., Марчевский И.К., Морева B.C. Моделирование обтекания тонкой пластинки с использованием модифицированной схемы метода вихревых элементов // Наука и образование. МГТУ им. Н.Э. Баумана. Электрон. журн. 2013. №9. C. 233-242. DOI: 10.7463/0913.0602362
11. Лаврентьев М.А., Шабат Б.В. Методы теории функций комплексного переменного. М.: Наука, 1965. 716 с.
12. Макарова М.Е. Поиск аналитических решений и исследование точности расчетных схем метода вихревых элементов в двумерных стационарных задачах обтекания профилей // Инженерный журнал: Наука и инновации: электронное научно-техническое издание.
2012. №4. Режим доступа: http://engjournal.ru/catalog/mathmodel/hidden/166.html (дата обращения 30.01.2015).
13. Морева B.C. Вычисление вихревого влияния в модифицированной схеме метода вихревых элементов // Инженерный журнал: Наука и инновации: электронное научно-техническое издание. 2012. №4. Режим доступа: http://engjournal.ru/search/author/ 268/page1.html (дата обращения 30.01.2015).
14. Морева B.C. Математическое моделирование обтекания профилей с использованием новых расчетных схем метода вихревых элементов: дис. .. . канд. физ.-мат. наук. М.,
2013. 130 с.
Science ¿¿Education
of the Bauman MSTU
Science and Education of the Bauman MSTU, 2015, no. 2, pp. 234-249.
DOI: 10.7463/0215.0756954
Received:
19.02.2015
Electronic journal
© Bauman Moscow State Technical University
ISSN 1994-0408
On the Accuracy of Numerical Schemes
for Flow Simulation Around Airfoils with Angle Point
Kuzmina K. S.1, Marchevsky I. K.1*, Moreva V. S.1
1 Bauman Moscow State Technical University, Russia
Keywords: convergence, integral equation, vortex element method, design model, singular solution
The paper considers two-dimensional problem of numerical flow simulation around an airfoil. The flow is assumed to be incompressible, so a meshless lagrangian vortex element method (VEM) can be used to simulate the flow. The aim of this research is to investigate the accuracy of the numerical schemes of VEM, which use different approaches to satisfy boundary conditions on the surface of airfoil in case of simulation of the flow around airfoils with angle point.
When constructing a numerical scheme, airfoil is replaced by polygon, which consists of panels. The airfoil is assumed to be rigid and immovable. The numerical schemes are based on two mathematically equivalent approaches to the boundary conditions satisfaction. The first ("classical") approach implies the equality to zero of the normal component of the flow velocity on the airfoil surface. It leads to the singular integral equation with Hilbert-type kernel, and the principal value of the corresponding integral should be understood in Cauchy sense. This equation should be solved numerically, but the mentioned features impose strong restrictions on discretization of airfoil, in particular, lengths of neighboring panels should not differ significantly. In the second case, the limit value of the tangential component of flow velocity from the airfoil side equates to zero. This approach leads to an integral equation of the second kind, which kernel is bounded in case of a smooth airfoil. This approach is free from limitations of the previous way. The numerical schemes based on the first and the second approaches are called NVEM and TVEM, correspondingly.
When solving engineering problems the airfoil sometimes should be split into panels with significantly different lengths. In order to investigate the accuracy of the mentioned NVEM and TVEM schemes the model problem of the flow simulation around the airfoil with angle point is considered. The asymptotic behavior of the solution, which is singular near the angle point, is known. Different algorithms of splitting airfoil into panels are used, and it is analyzed how it influences on the accuracy of the numerical schemes.
In case of the uniform split of the airfoil with angle point into panels both schemes NVEM and TVEM allow us to obtain correct solution with sufficiently high accuracy. If the split is non-uniform and panels, which are close to the angle point, have significantly different lengths, the NVEM scheme, in contrast to TVEM, doesn't allow us to obtain qualitatively correct solution.
References
1. Lifanov I.K. Metod singulyarnykh integral'nykh uravneniy i chislennyy eksperiment [Method of singular integral equations and numerical experiment]. Moscow, Yanus Publ., 1995. 521 p. (in Russian).
2. Andronov P.R., Guvernyuk S.V., Dynnikova G.Ya. Vikhrevye metody rascheta nestatsion-arnykh gidrodinamicheskikh nagruzok [Vortex methods of calculation of unsteady hydrody-namic loads]. Moscow, MSU Publ., 2006. 184 p. (in Russian).
3. Golovkin M.A., Golovkin V.A., Kalyavkin V.M. Voprosy vikhrevoi gidromekhaniki [Problems of vortex hydromechanics]. Moscow, Fizmatlit Publ., 2009. 264 p. (in Russian).
4. CottetG.-H., KoumoutsakosP.D. Vortex Methods: Theory andPractice. Cambridge University Press, 2008. 328 p.
5. Lewis R.I. Vortex Element Methods for Fluid Dynamic Analysis of Engineering Systems. Cambridge University Press, 2005. 592 p.
6. Kempka S.N., Glass M.W., Peery J.S., Strickland J.H. Accuracy Considerationsfor Implementing Velocity Boundary Conditions in Vorticity Formulations. SANDIA Report SAND96-0583 UC-700. Sandia National Laboratories, Albuquerque, N.M., March 1996. 52 p.
7. Moreva V.S., Marchevsky I.K. Vortex element method for 2D flow simulation with tangent velocity components on airfoil surface. ECCOMAS 2012 6th European Congress on Computational Methods in Applied Sciences and Engineering: Book of proceedings. Vienna, 2012. 14 p.
8. Kuzmina K.S., Marchevsky I.K. On Numerical Schemes in 2D Vortex Element Method for Flow Simulation around Moving and Deformable Airfoils. Advanced Problems in Mechanics (APM2014): Proc. ofthe XLII Summer School-Conference. St.-Petersburg, 2014, pp. 335-344.
9. Makarova M.E. Calculation of stationary attached flow past airfoil by stream of ideal incompressible medium. Vestnik MGTU im. N.E. Baumana. Ser. Estestvennye nauki = Herald of the BaumanMoscow State Technical University. Ser. Natural science, 2011, spec. iss. Prikladnaya matematika [Applied mathematics], pp. 63-74. (in Russian).
10. Makarova M.E., Marchevskii I.K., Moreva V.S. Flow simulation around a thin plate using a modified numerical scheme of the vortex element method. Nauka i obrazovanie MGTU im. N.E. Baumana = Science and Education of the Bauman MSTU, 2013, no. 9, pp. 233-242. DOI: 10.7463/0913.0602362 (in Russian).
11. Lavrent'ev M.A., Shabat B.V. Metody teorii funktsii kompleksnogo peremennogo [Methods of the theory of functions of a complex variable]. Moscow, Nauka Publ., 1965. 716 p. (in Russian).
12. Makarova M.E. Calculation of Flow past Airfoils of Simplest Topology Using the Modified Vortex-Element Method. Inzhenernyy zhurnal: nauka i innovatsii = Engineering Journal: Science and Innovation, 2012, no. 4. Available at: http://engjournal.ru/catalog/mathmodel/ hidden/166.html, accessed 30.01.2015. (in Russian).
13. Moreva V.S. Calculation of Vortex Effect in the Modified Numerical Scheme of Vortex-Element Method. Inzhenernyy zhurnal: nauka i innovatsii = Engineering Journal: Science and Innovation, 2012, no. 4. Available at: http://engjournal.ru/search/author/268/page1.html, accessed 30.01.2015. (in Russian).
14. Moreva V.S. Matematicheskoe modelirovanie obtekaniya profilei s ispol'zovaniem novykh raschetnykh skhem metoda vikhrevykh elementov. Kand. diss. [Mathematical modeling of flow of profiles using the new calculation schemes of vortex element method. Cand. diss.]. Moscow, 2013. 130 p. (in Russian).
