автореферат диссертации по информатике, вычислительной технике и управлению, 05.13.18, диссертация на тему:Математическое моделирование обтекания профилей с использованием новых расчетных схем метода вихревых элементов

кандидата физико-математических наук
Морева, Виктория Сергеевна
город
Москва
год
2013
специальность ВАК РФ
05.13.18
Диссертация по информатике, вычислительной технике и управлению на тему «Математическое моделирование обтекания профилей с использованием новых расчетных схем метода вихревых элементов»

Автореферат диссертации по теме "Математическое моделирование обтекания профилей с использованием новых расчетных схем метода вихревых элементов"

На правах рукописи

МОРЕВА Виктория Сергеевна

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

Специальность 05.13.18 — Математическое моделирование, численные методы и комплексы программ

АВТОРЕФЕРАТ диссертации на соискание ученой степени кандидата физико-математических наук

21 НОЯ 2013

005538784

Москва - 2013

005538784

Работа выполнена в федеральном государственном бюджетном образовательном учреждении высшего профессионального образования «Московский государственный технический университет имени Н.Э. Баумана»

Научный кандидат физико-математических наук, доцент

руководитель: Марчевский Илья Константинович

Официальные Галанин Михаил Павлович,

оппоненты: доктор физико-математических наук, профессор, федеральное государственное бюджетное учреждение науки «Институт прикладной математики имени М.В. Келдыша Российской академии наук», заведующий отделом

Дынникова Галина Яковлевна,

доктор физико-математических наук, старший научный сотрудник, Государственное учебно-научное учреждение «Научно-исследовательский институт механики Московского государственного университета имени М.В. Ломоносова», ведущий научный сотрудник

Ведущая федеральное государственное автономное образова-

организация: тельное учреждение высшего профессионального образования «Московский физико-технический институт (государственный университет)»

Защита состоится «10» декабря 2013 года в Ц час. 00 мин. на заседании диссертационного совета Д 212.141.15 при Московском государственном техническом университете им. Н.Э. Баумана по адресу: Москва, Рубцовская наб., д. 2/18, ауд. 1006 л.

С диссертацией можно ознакомиться в библиотеке Московского государственного технического университета им. Н.Э. Баумана.

Автореферат разослан «_»_2013 года.

Ученый секретарь диссертационного совета, кандидат технических наук, доцент

А.В. Аттетков

ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ

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

Бессеточные лагранжевы вихревые методы позволяют моделировать течения несжимаемой среды и решать широкий класс задач аэрогидродинамики и аэрогидроупругости. В основе вихревых методов лежит переход к завихренности как к первичной расчетной величине; непрерывное поле завихренности в области течения в расчете заменяется большим количеством вихревых элементов, движущихся по определенному закону в результате их парного взаимодействия. Значительный вклад в развитие вихревых методов внесли С.М. Белоцерковский и его последователи, А. Леонард, Г. Винкельманс, А. Чорин, Г. Котте, П. Комотсакос, К. Ка-мемото, Г.Я. Дынникова, Г.А. Щеглов, JI. Барба и др. Обоснованию вихревых методов посвящены работы И.К. Лифанова, A.B. Сетухи, Дж. Биля и А. Майды и др. Метод вязких вихревых доменов (П.Р. Андронов, C.B. Гувернюк, Г.Я. Дынникова) позволяет моделировать течения вязкой несжимаемой среды, однако требует значительных затрат вычислительных ресурсов, поэтому актуальной задачей является разработка алгоритмов и программ, позволяющих сократить затраты времени на проведение расчетов за счет использования параллельных вычислительных технологий и приближенных методов расчета парных взаимодействий вихревых элементов.

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

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

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

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

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

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

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

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

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

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

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

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

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

Научная новизна. В диссертации получены следующие новые научные результаты, которые выносятся на защиту.

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

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

3. Иерархия расчетных схем метода вихревых элементов, основанных на различных подходах к обеспечению выполнения граничного условия на профиле, с соответствующими расчетными формулами.

4. Параллельный программный комплекс РОЬАИА для расчета обтекания профилей методом вихревых элементов и вычисления действующих на них аэрогидродинамических нагрузок.

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

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

Разработан и зарегистрирован программный комплекс «РОЬАИА — Моделирование обтекания профилей и определение их аэродинамических нагрузок методом вихревых элементов» (свидетельство о государственной регистрации в Реестре программ для ЭВМ № 2013617462 от 14.08.2013 г.).

Апробация результатов работы. Результаты диссертационной работы апробированы на XIV, XV и XVI Международных симпозиумах «Методы дискретных особенностей в задачах математической физики» (г. Херсон, 2009, 2011 и 2013), V Всероссийской конференции «Необратимые процессы в природе и технике» (г. Москва, 2009), Междуна-

родной научной конференции «Параллельные вычислительные технологии» (г. Уфа, 2010 и г. Новосибирск, 2012), Всероссийской молодежной школе по параллельному программированию «Суперкомпьютерные технологии и высокопроизводительные вычисления в образовании, науке и промышленности» (г. Москва, 2010), V Международной конференции «International Conference on Vortex Flows and Vortex Models» (r. Ka-зерта, Италия, 2010), X Всероссийском съезде по фундаментальным проблемам теоретической и прикладной механики (г. Нижний Новгород, 2011), XXXVIII Международной молодежной научной конференции «Гагаринские чтения» (г. Москва, 2012), Международной конференции «International Conference on Computational Fluid Dynamics» (г. Париж, 2012), VI Международном конгрессе «Computational Methods in Applied Sciences and Engineering (ECCOMAS-2012)» (г. Вена, 2012), Международной летней школе-конференции «Advanced Problems in Mechanics» (г. Санкт-Петербург, 2013). Результаты исследований обсуждались также на Международном авиационно-космическом научно-гуманитарном семинаре им. С. М. Белоцерковского (ЦАГИ им. проф. Н.Е. Жуковского, ВВИА им. Н.Е. Жуковского и Ю. А. Гагарина, 2010, 2012).

Диссертация является составной частью фундаментальных исследований, проводимых в рамках гранта Президента РФ для государственной поддержки ведущих научных школ (грант НШ-255.2012.8) и гранта Президента РФ для государственной поддержки молодых российских ученых — кандидатов наук (грант МК-6482.2012.8).

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

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

Структура и объем работы. Диссертация состоит из введения, трех глав, заключения и списка литературы. Диссертационная работа изложена на 130 страницах, содержит 42 иллюстрации и 13 таблиц. Библиография включает 92 наименования.

СОДЕРЖАНИЕ РАБОТЫ

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

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

dV Vp

V • V = 0, -—- + (V • V)F = uAV---, reS.

dt p

Здесь V = V{r, t) — скорость среды в точке г(х,у) в момент времени i; р = p(r, t) — давление; р = const — плотность среды; и — коэффициент кинематической вязкости; 5 — область течения.

На неподвижном обтекаемом профиле К задается граничное условие прилипания, а на бесконечности — условия затухания возмущений:

re К: V(r,t)= О, |г| —► оо: V{r,t)^Voa, p{r,t)-+p0о.

Для решения поставленной задачи используется метод вихревых элементов (МВЭ), в котором первичной расчетной величиной является завихренность ÎÎ = V х V, по известному распределению которой с помощью закона Био — Савара может быть восстановлено поле скоростей V. Распределение завихренности моделируется набором из N вихревых элементов — вихревых нитей, перпендикулярных области течения, характеризуемых положениями в потоке Г; и циркуляциями Г,-, i = 1.....N.

Скорость среды V в точке г в фиксированный момент времени вычисляется по формуле

ViD-V. + gQir-r,)*,

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

Щп) = vx + V(n) + Wird, W(n) = где W — диффузионная скорость пропорциональная вязкости; П = fî • k.

Для вычисления действующих на профиль нагрузок используются аналог интеграла Коши — Лагранжа, позволяющий найти распределение давления по профилю, или формула для расчета главного вектора и главного момента аэродинамических сил, разработанные Г.Я. Дынниковой (2000).

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

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

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

Традиционно при решении задачи Коши для указанной выше системы обыкновенных дифференциальных уравнений используется явный метод Эйлера, поэтому точность можно повысить, если использовать, например, метод Рунге — Кутты второго порядка точности. Несмотря на имеющиеся упоминания о таком подходе (Г.А. Павловец, A.C. Петров, 1974; С.М. Белоцерковский и др., 1997), систематических исследований в этом направлении до последнего времени не проводилось.

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

Эффективность применения метода Рунге — Кутты показана на примере модельной задачи о расчете диффузии вихря в вязкой жидкости. На Рис. 1 представлены результаты вычисления циркуляции поля скоростей вдоль окружности фиксированного радиуса, полученные с применением методов Эйлера и Рунге — Кутты 2-го порядка, и их сравнение с аналитическим решением.

8000 12000 16000 8000 12000 16000 Рис. 1. Циркуляция поля скоростей при эволюции вихря в идеальной (а) и вязкой (б) средах

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

Во второй главе рассмотрены два возможных подхода МВЭ к построению математической модели генерации завихренности на профиле. В основе вихревых методов лежит идея Н. Е. Жуковского о замене профиля вихревым слоем определенной интенсивности. При математическом моделировании течений вязкой среды методом вихревых элементов данный вихревой слой не является присоединенным и полностью сходит в поток. Интенсивность вихревого слоя 7 = 7 • k на профиле может быть найдена из граничного условия на нем.

Предельное значение скорости среды со стороны профиля равно

-'■* I - * ») - ™ - ^

где п(г) — вектор единичной внешней нормали к профилю; т(г) — единичный вектор касательной в точке г е К, при этом я(г) х т(г) = к.

Для нахождения неизвестной интенсивности вихревого слоя 7(г) используется условие равенства нулю вектора скорости V_(r) на профиле К — условие прилипания. Если вся генерируемая завихренность становится свободной, то это условие обеспечивается равенством нулю компоненты вектора V_(r) на профиле: либо нормальной (подход НМВЭ), либо касательной (подход КМВЭ, S.N. Kempka et al., 1996). Использование подхода НМВЭ приводит к сингулярному интегральному уравнению с неограниченным ядром, интеграл в котором понимается в смысле главного значения по Коши:

V.(r) ■ п{г) = 0 <=> £ [k dlro = -и(г) ■

Подход КМВЭ приводит к интегральному уравнению Фредгольма второго рода с ограниченным (для гладкого профиля) ядром:

У_(г) • т(г) - 0 £ [*Х^,;:УГ)7(гоЖ0 - ^ = -Пг) ■ V«,.

Для выделения единственного решения уравнений используется аналог теоремы Томсона для вязкой среды (Г.Я. Дынникова, 2003), согласно которой суммарная завихренность в области течения сохраняется. В случае нулевой суммарной завихренности в начальный момент времени получаем: = 0. к

Для решения задачи методом вихревых элементов разработаны восемь расчетных схем, которые получаются различными комбинациями следующих подходов (Таблица 1):

• используется подход НМВЭ или КМВЭ;

• граничное условие (ГУ) выполняется в контрольных точках (точках коллокации) или в среднем на панелях профиля;

• завихренность на каждой панели стянута в вихревой элемент (ВЭ) или распределена по панели.

Таблица 1

НМВЭ (АО / КМВЭ (Г) подходы ГУ в точках коллокации (,coll) ГУ в среднем на панели (aver)

Вихревые элементы (vort) kt'coll 1 q-coll Jvvort l vort kfaver t Taver •"'vort / vort

Вихревой слой (layer) kfcoll I fcoll •/4laver ' layer kfaver / fever J layer ' latter

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

Ыр ( (\ \\ Ыр £ (А„ - Д ( ) ) Чі + = -Ві (і = 1,.. .Мр), =

7=1 ^ V // /=1

Величины, входящие в данную СЛАУ, имеют следующий смысл:

Ац — коэффициент матрицы, учитывающий влияние завихренности с у-й панели на г'-ю панель: для подхода НМВЭ необходимо вычислять все коэффициенты Ац, для подхода КМВЭ Ац = 0; Д — диагональные элементы матрицы: Д = 1 для схем Тщ^ и

Д = 1/Д/і для схем Тис°г!1 и где Д/,- — длина г-й панели, Д = 0

для всех схем НМВЭ; 5ц — символ Кронекера, 5ц = 1 при г = у и <5(/- = 0 при і ф у; Ьі — коэффициенты при регуляризирующей переменной: для схем НМВЭ Ьі = Мі, для схем КМВЭ Ьі = 1;

М] — коэффициенты уравнения, аппроксимирующего дополнительное условие: для схем, в которых завихренность стягивается в вихревые элементы М) = 1; для схем с распределенной по панелям завихренностью М) выражается через длины панелей;

В,- — проекции суммарной скорости набегающего потока и влияния вихревого следа на нормаль к г-й панели — для схем подхода НМВЭ и на касательную к ¿-й панели — для подхода КМВЭ.

Неизвестными величинами в системе являются:

qj — циркуляции вихревых элементов или значения кусочно-постоянной интенсивности вихревого слоя; — регуляризирующая переменная (И.К. Лифанов, 1995).

В диссертации для всех рассмотренных расчетных схем получены формулы для вычисления вышеперечисленных коэффициентов и правых частей СЛАУ, которые могут быть включены в соответствующие компьютерные программы.

На примере модельных задач о расчете стационарного обтекания профилей простых форм с известным аналитическим решением выполнено сравнение результатов расчетов, проведенных с использованием рассмотренных расчетных схем, с известным аналитическим решением. Величина погрешности указывает на максимальную ошибку в определении циркуляции вихревых элементов, которые генерируются на профиле и затем сходят в поток. Погрешности схем, основанных на подходах НМВЭ и КМВЭ, вычисляются по формулам

||ДГ||лг = ¿тах[(|7; - т?|)(Д/«-1 + А«]. Д'о = Д/„;

||ДГ||Т = тах[(|7г-7,0|)А/1],

где 7,- — найденное в расчетах среднее значение интенсивности вихревого слоя на ¿-й панели, 7° — среднее значение точной величины интенсивности на рассматриваемом элементе, А1— длина панели.

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

Таблица 2

||ДГ||лг / ||ДГ||Т Эллиптический профиль Профиль Жуковского

д[call j rrcolt Jvvort 1 vort 0,0058 / 0,0055 4,0815/0,0108

Kfcoll j fcoll J laifer ' lauer 0,0013 / 0,0071 0,0398 / 0,0163

Kfaver ! *j~aver Jvvort 1 vor t 0,0182 / 0,0009 0,0544 / 0,0158

Kfaver 1 rr~üver lauer < 1lauer 0,0008 / 0,0006 0,0011 /0,0004

На примере расчета обтекания профилей с тремя угловыми точками показано, что использование подхода КМВЭ приводит к результатам, хорошо согласующимся с асимптотическим решением (И.К. Лифа-нов, 1995). На Рис. 2 показан пример такого профиля при равномерном и неравномерном его разбиении на панели вблизи угловой точки.

■ а, в '

Рис. 2. Профиль с угловыми точками (ф = 7п/6, а = 7г/6), равномерное разбиение (а), неравномерное разбиение (б)

При одинаковых длинах панелей, примыкающих к угловой точке (равномерное разбиение), решения, получаемые при использовании всех расчетных схем, хорошо согласуются с асимптотическим. При различающихся длинах панелей (неравномерное разбиение) результаты расчета для схем М™уеегг и %с°'} приведены на Рис. 3. В отличие от схемы Л/"™", решение, полученное с использованием схемы близко к асимптотическому.

7.0 5.0 2.0

Л5Й

3.5 3.0 2.5

•7-aver 'layer

0.001 0.002 0.003

0.001 0.002 0.003

Рис. 3. Результаты расчета средней интенсивности вихревого слоя на панелях при неравномерном разбиении профиля (точки — результат расчета, сплошная линия — асимптотическое решение)

На примере решения задачи о математическом моделировании стационарного обтекания кругового профиля с изолированным вихрем вблизи него показано, что наибольшая точность вычислений получается при использовании расчетной схемы Т^". В Таблице 3 показана погрешность вычислений при расстоянии от вихря до профиля с1\ = 0,01 и ¿2 = 0,05.

Таблица 3

||ДГ|М|ДГ||Г di =0,01 ¿2 = 0,05

Kfcoll / H~COll 0,1251 / 0,1754 0,0136 / 0,0046

Kfcoll / '■rcoll layer / lauer 0,0326 / 0,1758 0,0250 / 0,0051

kfaver / ¿-raver ' vorl t vort 0,2178 / 0,0017 0,0642 / 0,0007

Kfaver j q-aver Jylai/er ' lauer 0,1842 / 0,0016 0,0158 / 0,0006

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

В третьей главе диссертации описан разработанный программный комплекс РОЬАИА, позволяющий решать серии однотипных нестационарных задач по математическому моделированию процесса обтекания профилей потоком вязкой несжимаемой среды. В его основу положены алгоритм МВЭ, описанный в главе 1, с расчетными схемами, представленными в главе 2. В каждой из задач моделируется обтекание профиля под фиксированным углом атаки. Использование параллельных алгоритмов позволяет автоматизировать процесс решения поставленных задач: по мере завершения одних расчетов и освобождения соответствующих вычислительных ядер запускаются процедуры решения следующих задач. Предусмотрена возможность распараллеливания наиболее трудоемких операций алгоритма. Программный комплекс способен работать на различных типах вычислительных систем — от персональных компьютеров до высокопроизводительных кластеров.

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

Рис. 4. Схема распараллеливания серии задач

Этап 1. На этапе запуска расчета из Ыргос процессоров выделяется головной процессор, на который загружается список углов атаки {а;}^"4 с указанием чисел Пр)ос, необходимых для моделирования обтекания при соответствующем угле атаки. Глобальный головной процессор контролирует решение серии задач и осуществляет пересылки необходимых параметров между другими процессорами.

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

Этап 3. На данном этапе непосредственно выполняются шаги расчета обтекания профиля с помощью МВЭ, операции которого отображены на Рис. 5. Расчет повторяется, пока не закончится квант времени или выполнится критерий прекращения счета.

Рис. 5. Схема выполнения шага расчета методом вихревых элементов

Для оценки эффективности использования параллельных алгоритмов вычислено ускорение при решении серии из 91 задачи по расчету обтекания профиля крыла ЦАГИ Р2-18 на вычислительной машине МВС-100К МСЦ РАН (Таблица 4).

Таблица 4

Число ядер 1 2 4 8 16 32 48 64

Ускорение (1 ядро/расчет) 1,00 1,93 3,76 7,21 14,04 27,22 39,43 47,33

Ускорение (2 ядра/расчет) — 1,92 3,74 7,37 14,30 28,28 41,28 53,79

Также для ускорения расчетов в программном комплексе РОЬАИА реализован быстрый метод приближенного расчета парных взаимодействий вихревых элементов, описанный Г.Я. Дынниковой (2009). Данный метод основан на построении дерева подобластей, содержащих вихревые элементы. Получена уточненная теоретическая оценка вычислительной трудоемкости алгоритма быстрого метода:

0 =

896тг У Я4

1п(р- 1) -

р- 1

+ 1

+ 4АГ +

24ЛГ2

р = 9 ■ (л/2)*"4,

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

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

На Рис. 6-7 приведены зависимости аэродинамических коэффициентов кругового профиля от номера шага расчета, полученные с использованием схем (с осреднением и без него) и на начальном этапе развития вихревого следа. Использование схемы приводит к существенно меньшим осцилляциям характеристик от шага к шагу.

200 400 600 800

(с осреднением по методу скользящей средней)

Номер шага по времени

0 200 400 600 800

Номер шага по времени

0 200 400 600 800

Рис. 6. Нестационарный коэффициент лобового сопротивления круга

(с осреднением по методу скользящей средней)

Номер шага по времени

200 400 600 800

0 200 400 600 800

Номер шага но времени

О 200 400 600 800

Рис. 7. Нестационарный коэффициент подъемной силы круга

Зависимости коэффициентов лобового сопротивления и подъемной силы кругового профиля от времени показаны на Рис. 8.

1.5 1.0 0.5 0.0

1.0 о.о -1.0 -2.0,

ело

10 20 30 40 50 60

10 20 30 40 50 60

Рис. 8. Нестационарные аэродинамические коэффициенты круга

Среднее значение коэффициента лобового сопротивления при t > 30 составляет Сха « 1,19, амплитуда колебаний коэффициента подъемной силы Суа и 1,5, безразмерная частота схода вихрей БЬ » 0,22. Полученные результаты хорошо согласуются с данными экспериментов.

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

-pj.ee)

30 60 90 120 150 180 -2

Рис. 9. Зависимости стационарных коэффициентов лобового сопротивления и подъемной силы полукруглого профиля от угла атаки (точки — результат расчета, сплошная линия — экспериментальные данные (З.П. Случановская, 1973))

Для задачи Блазиуса о моделировании обтекания тонкой пластинки вид близкого к установившемуся вихревого следа представлен на Рис. 10.

ц-]-1

Рис. 10. Вихревой след вблизи пластинки

На Рис. 11 в безразмерных величинах представлены распределения продольной и поперечной компонент скорости в сечении х = 0,35. 1.0

Рис. 11. Графики распределения продольной и поперечной компонент скорости в сечении х = 0,35 пограничного слоя на пластинке (точки — расчет, сплошная линия — решение Блазиуса)

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

ОСНОВНЫЕ РЕЗУЛЬТАТЫ ДИССЕРТАЦИОННОЙ РАБОТЫ

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

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

3. Разработан параллельный программный комплекс РОЬА1?А, позволяющий решать серии задач и проводить отдельные расчеты на многопроцессорных вычислительных системах по математическому моделированию нестационарного обтекания профилей и определению их аэродинамических характеристик.

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

ОСНОВНЫЕ РЕЗУЛЬТАТЫ ДИССЕРТАЦИИ ОТРАЖЕНЫ В РАБОТАХ

1. Марчевский И.К., Морева B.C. Уточненный метод вихревых элементов для моделирования обтекания профиля потоком // Необратимые процессы в природе и технике: труды V Всероссийской конференции. М., 2009. С. 97-101.

2. Морева B.C. Применение метода второго порядка точности при моделировании течений методом вихревых элементов // Методы дискретных особенностей в задачах математической физики: труды XIV Международного симпозиума. Харьков-Херсон, 2009. С. 133-136.

3. Марчевский И.К., Морева B.C. Высокопроизводительный программный комплекс POLARA для определения аэродинамических характеристик профилей // Параллельные вычислительные технологии: труды Международной научной конференции. Уфа, 2010. С. 533-538.

4. Марчевский И.К., Морева B.C. Численное моделирование обтекания системы профилей методом вихревых элементов // Вестн. Моск. гос. техн. ун-та им. Н.Э. Баумана. Сер.: Естественные науки. 2010. № 1. С. 12-20.

5. Moreva V.S., Marchevsky I.К. High-efficiency POLARA program for airfoil aerodynamic characteristics calculation using vortex elements method // The 5th International Conference on Vortex Flows and Vortex Models: book of proceedings. Caserta (Italy), 2010. P. 1-6.

6. Марчевский И.К., Морева B.C. Параллельный программный комплекс POLARA для исследования расчетных схем метода вихревых элементов // Методы дискретных особенностей в задачах математической физики: труды XV Международного симпозиума. Харьков-Херсон, 2011. С. 280-283.

7. Исследование аэроупругих колебаний провода, вызываемых отрывным вихревым обтеканием / B.C. Морева [и др.] // Вестник Нижегородского университета им. Н.И. Лобачевского. 2011. Ч. 2. № 4. С. 157-159.

8. Морева B.C. Способы ускорения вычислений в методе вихревых элементов // Вестн. Моск. гос. техн. ун-та им. Н.Э. Баумана. Сер.: Естественные науки. Спец. выпуск «Прикладная математика». 2011. С. 83-95.

9. Марчевский И.К., Морева B.C. Параллельный программный комплекс POLARA для моделирования обтекания профилей и исследования расчетных схем метода вихревых элементов // Параллельные вычислительные технологии: труды Международной научной конференции. Новосибирск, 2012. С. 236-247.

10. Марчевский И.К., Морева B.C. Численное моделирование в задачах аэрогидродинамики с использованием метода вихревых элементов // CAD-CAM-CAE Observer. 2012. № 2. С. 84-91.

11. Marchevsky I.К., Moreva V.S. On modified numerical schemes in vortex element method for 2D flow simulation around airfoils // World Academy of Science, Engineering and Technology. 2012. V. 68. P. 1594-1600.

12. Moreva V.S., Marchevsky I.K. Vortex element method for 2D flow simulation with tangent velocity components on airfoil surface // 6th European Congress on Computational Methods in Applied Sciences and Engineering: book of proceedings. Vienna, 2012. 14 p.

13. Морева B.C. Вычисления вихревого влияния в модифицированной численной схеме метода вихревых элементов // Вестн. Моск. гос. техн. ун-та им. Н.Э. Баумана. Сер.: Естественные науки. 2012. Спец. выпуск № 2 «Математическое моделирование в технике». С. 137-144.

14. Учебно-экспериментальный вычислительный кластер. 4.2. Примеры решениязадач / B.C. Морева [и др.] // Вестн. Моск. гос. техн. ун-та им. Н.Э. Баумана. Сер.: Естественные науки. 2012. № 4. С. 82-102.

15. Макарова М.Е., Марчевский И.К., Морева B.C. О различных подходах к моделированию обтекания профиля методом вихревых элементов // Методы дискретных особенностей в задачах математической физики: труды XVI Международного симпозиума. Харьков-Херсон, 2013. С. 246-249.

16. Moreva V.S. On simulation of the flow around an airfoil using different numerical schemes of vortex element method // Advanced Problems in Mechanics: proceedings of the XLI Summer School-Conference. St.Petersburg, 2013. P. 387-396.

17. Свидетельство о государственной регистрации программы для ЭВМ № 2013617462. POLARA — Моделирование обтекания профилей и определение их аэродинамических нагрузок методом вихревых элементов / B.C. Морева, И.К. Марчевский, Г.А. Щеглов. Зарегистрировано в Реестре программ для ЭВМ 14.08.2013.

18. Макарова М.Е., Марчевский И.К., Морева B.C. Моделирование обтекания тонкой пластинки с использованием модифицированной схемы метода вихревых элементов // Наука и образование: электронное научно-техническое издание. 2013. № 9. DOI: 10.7463/0913.0602362.

Режим доступа: http://technomag.bmstu.ru/doc/602362.html.

Подписано к печати 1.11.13. Заказ № 702 Объем 1,0 печ.л. Тираж 100 экз. Типография МГТУ им. Н.Э. Баумана 105005, Москва, 2-я Бауманская ул., д.5,стр.1 (499) 263-62-01

Текст работы Морева, Виктория Сергеевна, диссертация по теме Математическое моделирование, численные методы и комплексы программ

МОСКОВСКИЙ ГОСУДАРСТВЕННЫЙ ТЕХНИЧЕСКИЙ УНИВЕРСИТЕТ имени Н.Э. Баумана

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

05.13.18 — Математическое моделирование, численные

0420136495

На правах рукописи

Морева Виктория Сергеевна

методы и комплексы программ

Диссертация на соискание ученой степени кандидата физико-математических наук

Научный руководитель: к.ф.-м.н., доцент И.К. Марчевский

Москва - 2013

Оглавление

Стр.

Введение............................. 5

Глава 1. Метод вихревых элементов для моделирования течений вязкой несжимаемой среды............. 17

1.1. Математическая постановка задачи............. 17

1.2. Восстановление поля скоростей по известному распределению завихренности и расчет аэродинамических нагрузок 21

1.3. Моделирование распределения завихренности с помощью вихревых элементов...................... 25

1.4. Расчет движения вихревых элементов и вычисление аэродинамических нагрузок.................... 26

1.5. Метод второго порядка точности.............. 30

1.6. Верификация метода Рунге — Кутты второго порядка точности .............................. 33

1.6.1. Математическое моделирование диффузии круглого вихря........................... 34

1.6.2. Математическое моделирование обтекания кругового профиля при различных значениях числа Рей-нольдса.......................... 36

1.7. Результаты главы 1...................... 37

Глава 2. Расчетные схемы метода вихревых элементов ... 39

2.1. Подходы к моделированию генерации завихренности на обтекаемом профиле...................... 39

2.1.1. Равенство нулю нормальной компоненты скорости на профиле......................... 41

2.1.2. Равенство нулю касательной компоненты скорости

на профиле........................ 45

2.2. Классическая расчетная схема на обтекаемом профиле . . 47

Стр.

2.3. Новые расчетные схемы метода вихревых элементов ... 50

2.4. Коэффициенты расчетных схем метода вихревых элементов 53

2.4.1. Расчетные схемы и ..........................56

2.4.2. Расчетные схемы и Г™ге{........................57

2.4.3. Расчетные схемы Мс^ег и ........................59

2.4.4. Расчетные схемы и ........................61

2.5. Определение циркуляций вихревых элементов, сходящих

в поток..........................................................64

2.6. Верификация расчетных схем метода вихревых элементов 65

2.6.1. Математическое моделирование стационарного обтекания эллиптического профиля и профиля Жуковского ......................................................67

2.6.2. Математическое моделирование стационарного обтекания профиля с тремя угловыми точками............71

2.6.3. Математическое моделирование обтекания кругового профиля с изолированным вихрем..................75

2.7. Результаты главы 2............................................77

Глава 3. Программный комплекс РОЬАЕА....................79

3.1. Структура программного комплекса РОЬАИА..............79

3.2. Оценка эффективности использования параллельных алгоритмов ........................................................87

3.2.1. Оценка эффективности внутреннего распараллеливания ......................................................88

3.2.2. Оценка эффективности внешнего распараллеливания 90

3.2.3. Исследование возможности применения различных моделей параллельного программирования............93

3.3. Быстрый метод расчета парных взаимодействий вихревых

элементов........................................................94

3.3.1. Описание метода и оценка его эффективности ... 94

Стр.

3.3.2. Оценка оптимальной глубины дерева в быстром методе 99

3.3.3. Параллельная реализация быстрого метода..... 103

3.4. Использование программного комплекса POLARA для проведения методических исследований параметров расчетных схем метода вихревых элементов.......... 105

3.5. Верификация программного комплекса POLARA..... 108

3.5.1. Моделирование обтекания профиля при использовании расчетных схем //¡¡Ц и ........... 108

3.5.2. Моделирование обтекания кругового профиля при

Re = 103......................... 110

3.5.3. Задача Блазиуса..................... 113

3.5.4. Моделирование обтекания профиля в форме полукруга117

3.5.5.Моделирование обтекания квадратного профиля . . 118

3.6. Результаты главы 3...................... 119

Основные результаты и выводы................120

Литература...........................121

Введение

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

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

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

Первые систематические исследования в области математического моделирования обтекания профилей и вычисления действующих на них нагрузок были проведены в начале XX века Н.Е. Жуковским [17]. С использованием аппарата теории функций комплексного переменного им было получено решение стационарной задачи о расчете безотрывного обтекания профилей простых форм (для которых удается построить конформное отображение внешности профиля на внешность круга) неогра-

ниченным потоком идеальной несжимаемой среды. Именно Н.Е. Жуковским установлена возможность моделирования профиля в потоке при помощи вихревого слоя, расположенного на профиле. Им же в 1904 г. была получена известная формула, в соответствии с которой подъемная сила крыла пропорциональна циркуляции поля скоростей вдоль профиля. Аналитически величину циркуляции удается вычислить для крыловых профилей в предположении, что поток сходит с острой кромки; данное предположение положено в основу гипотезы Чаплыгина — Жуковского, которая впоследствии была обобщена и в настоящее время более известна как гипотеза Чаплыгина — Жуковского — Кутты. В отечественной литературе исследования в данном направлении получили значительное развитие в трудах М. В. Келдыша, Л. И. Седова, М. А. Лаврентьева [23, 28, 63] и др.

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

ревые методы, позволяющие моделировать течения несжимаемой среды, решать широкий класс задач аэрогидродинамики и аэрогидроупругости с приемлемой точностью и требующие меньших затрат вычислительных ресурсов. В основе вихревых методов лежит переход к завихренности как к первичной расчетной величине; непрерывное поле завихренности в области течения в расчете заменяется большим количеством вихревых элементов, движущихся по определенному закону в результате их парного взаимодействия. Вычислительные ресурсы таким образом сосредоточены в той части области течения, в которой завихренность отлична от нуля; эта область обычно имеет сравнительно небольшие размеры. Значительный вклад в развитие вихревых методов внесли научные школы С. М. Белоцерковского [66], А. Леонарда [78], Г. Винкельман-са [91], Г. Котте [76], К. Камемото [90], Л. Барбы [72] и др. Среди работ по обоснованию вихревых методов можно выделить исследования И. К. Лифанова [5], А. В. Сетухи [24], Дж. Биля и А. Майды [74]. Достаточно полный обзор вихревых методов имеется в работах Т. Сарпкайи и А. Леонарда [60, 78].

Метод дискретных вихрей [30, 57, 66] позволяет моделировать отрывное обтекание профиля, имеющего угловые точки, потоком идеальной несжимаемой среды. Вычислительная трудоемкость метода достаточно низкая, что способствовало широкому применению метода и позволило решить многие практически важные задачи, используя весьма скромную по современным меркам вычислительную технику 70-х-80-х годов XX века. Метод дискретных вихрей оказался исключительно эффективным при моделировании обтекания элементов конструкций летательных аппаратов. К недостаткам классического метода дискретных вихрей следует отнести невозможность моделирования обтекания гладких профилей, точка отрыва потока с которых неизвестна, а также отрыва потока с гладких участков профилей, например, крыловых, при больших углах атаки. Известны попытки преодоления этих ограничений метода

дискретных вихрей: точки отрыва потока предлагалось находить путем построения приближенного решения уравнений пограничного слоя [47] или используя некоторые вариационные принципы [9]. Данные подходы в силу их сложности и неуниверсальности широкого распространения не получили. Также важно отметить, что в рамках метода дискретных вихрей можно моделировать только течения идеальной среды, тогда как на практике возникает необходимость расчета течений вязкой среды, характеризуемых умеренными значениями числа Рейнольдса, когда вязкие эффекты оказывают значительное влияние на структуру течения и аэродинамические нагрузки.

Попытки преодоления отмеченных недостатков метода дискретных вихрей при моделировании двумерных (плоских и осесимметричных) течений привели к разработке модификаций вихревых методов, позволяющих приближенно моделировать отрывное обтекание гладких профилей — этому вопросу посвящены работы А. Чорина [75], а также Г. А. Пав-ловца и A.C. Петрова [54, 55]. Вычислительные ресурсы ЭВМ 70-х годов, однако, не позволили развить эти подходы. Эти идеи, в частности, были положены в основу модификации метода дискретных вихрей, разработанной Г. А. Щегловым и И. К. Марчевским [37, 46].

Также были созданы вихревые лагранжевы методы расчета вязких течений, описываемых уравнениями Навье — Стокса. Известно несколько подходов к учету влияния вязкости, наибольшее распространение получили метод случайных блужданий, предложенный А. Чориным [75] и развитый другими исследователями (достаточно подробный обзор работ в этом направлении содержится в диссертации [88]) и метод диффузионный скорости, предложенный И. Огами и Т. Акаматсу [84] и значительно развитый Г.Я. Дынниковой [13, 15, 16]. В работах [1, 12] описан разработанный сотрудниками Института механики МГУ им. М. В. Ломоносова метод вязких вихревых доменов, а также эффективные способы определения нагрузок, действующих на обтекаемые профили, являющи-

еся обобщениями интеграла Коши — Лагранжа на случай вихревого движения вязкой несжимаемой среды.

В отличие от классического метода дискретных вихрей, реализация метода вязких вихревых доменов требует существенно больших вычислительных ресурсов. Количество вихревых элементов, моделирующих вихревой след за профилем, на практике может достигать десятков и сотен тысяч, поэтому непосредственный расчет их парных взаимодействий, даже на современных ЭВМ, приводит к недопустимо большим затратам машинного времени: расчет нестационарного обтекания профиля до выхода на установившийся режим на персональных ЭВМ может продолжаться несколько суток и даже недель. Поэтому актуальным вопросом является разработка алгоритмов, позволяющих сократить затраты времени на проведение расчетов. Значительный эффект дает использование параллельных вычислительных технологий, поскольку наиболее затратная часть алгоритма — расчет парных взаимодействий вихревых элементов — хорошо распараллеливается. Однако в известных работах по распараллеливанию метода вихревых элементов авторами не уделялось внимание распараллеливанию других этапов алгоритма. В частности, в работе [71] получено практически линейное ускорение на этапе расчета парных взаимодействий, однако если задействовать в расчете более, чем 4-8 вычислительных ядер, то эффективность распараллеливания при решении задачи в целом резко снижается.

Иной подход к ускорению вычислений предполагает использование приближенных методов расчета парных взаимодействий вихревых элементов, вычислительная трудоемкость которых составляет 0(ЛПо§Л/") против 0(Ы2) при «прямом» расчете. К ним относятся методы мульти-польных разложений, например, аналогичные методу Барнса — Хата быстрого решения «задачи N тел» [73], методы скелетонных разложений [67], а также некоторые другие подходы, обзор которых можно найти в диссертации [89].

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

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

Использование эффективных программных реализаций позволяет при сохранении приемлемыми временных затрат существенно увеличивать количество вихревых элементов, моделирующих обтекаемый профиль и вихревой след. Однако, как показывают систематические расчеты, с ростом количества вихревых элементов точность моделирования возрастает далеко не всегда. Это особенно актуально при расчете обтекания профилей, имеющих острую кромку: вблизи кромки интенсивность вихревого слоя на профиле определяется со значительной погрешностью, и эта погрешность возрастает с увеличением количества вихревых элементов, моделирующих профиль. О данной проблеме упоминается, в частности, в [30], там же показано, что при выполнении определенных условий к гладкости профиля и поля скоростей среды численное решение задачи, получаемое методом, который традиционно используется в классическом методе дискретных вихрей, сходится к точному решению

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