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

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

Плёнкин Андрей Валерьевич

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

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

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

005546476

2 7 НАР гон

Москва - 2013

005546476

Работа выполнена на кафедре вычислительной механики механико-математического факультета Федерального государственного образовательного учреждения высшего профессионального образования "Московский государственный университет имени М.В.Ломоносова".

Научный руководитель:

доктор физико-математических наук,

Афендиков Андрей Леонидович

Официальные оппоненты:

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

Георгиевский Павел Юрьевич, кандидат физико-математических наук, Федеральное государственное бюджетное образовательное учреждение высшего профессионального образования "Московский государственный университет имени М.В.Ломоносова", Научно-исследовательский институт механики МГУ, ведущий научный сотрудник

Ведущая организация:

Федеральное государственное бюджетное учреждение науки Институт вычислительной математики Российской академии наук (ИВМРАН).

Защита состоится " L Ч " 2014 г. в 7"Л' О О час. на

заседании Диссертационного совета Д 002.024.03, созданного на базе Федерального государственного бюджетного учреждения науки Институт прикладной математики им. М.В. Келдыша Российской академии наук (ИПМРАН) по адресу: 125047, Москва, Миусская пл., 4.

С диссертацией можно ознакомиться в библиотеке ИПМ им. М.В. Келдыша

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

РАН.

Автореферат разослан

к

Ученый секретарь диссертационного совета Д 002.024.03, доктор физико-математических наук Змитренко Н.В.

Общая характеристика работы

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

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

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

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

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

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

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

Исследования, вошедшие в диссертацию, были частично поддержаны грантами РФФИ 08-01-00454-а, 11-01-00390-а и программой 3 ОМН РАН.

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

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

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

В первую очередь следует упомянуть монографию H.H. Яненко и Е.В. Ворожцова, а также развивающие эти работы исследования С Б. Базарова. Кроме того имеется значительное количество работ, посвященных решению задачи выделения краев,

4

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

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

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

Положения, выносимые на защиту:

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

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

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

4) Проведено численное моделирование и выделение разрывов в ряде прикладных задач газовой динамики. Рассмотрены сверхзвуковое обтекание тел под углом атаки и распространение ударных волн в каналах при наличии импульсного энерговыделения. Расчеты проводились для идеальных и вязких течений при больших числах Рейнольдса. Показано, что «разрывы» выделяемые в найденных численно решениях уравнений Навье-Стокса при больших числах Рейнольдса вдали от тела близки к разрывам в решениях уравнений Эйлера. В решениях уравнений Навье-Стокса алгоритмически выделяются погранслой и вихревые структуры свойственные только вязким течениям.

Степень достоверности и апробация результатов. Основные результаты диссертации докладывались на следующих научно-исследовательских семинарах и конференциях: научно-исследовательские семинары кафедры вычислительной механики механико-математического факультета Московского государственного университета им. М.В.Ломоносова (МГУ им. М.В.Ломоносова 2009-2011 гг.); научно-исследовательские семинары Института прикладной математики им. М.В. Келдыша РАН (ИВМ РАН, 2011 г.); научно-исследовательский семинар Института вычислительной математики РАН (ИВМ РАН, 2013 г.); XVI Международная конференция студентов, аспирантов и молодых учёных «Ломоносов» (г. Москва 2009 г.); научная конференция «Ломоносовские чтения» (г. Москва 2009 г.); XVIII (2010 г.), XIX (2012 г.) всероссийские конференции «Теоретические основы и конструирование численных алгоритмов решения задач математической физики», посвященная памяти К.И. Бабенко (Абрау-Дюрсо, Новороссийск, Россия); XIV (2011 г.), XV (2013 г.) всероссийские молодежные конференции-школы с международным участием «Современные проблемы математического моделирования» (Абрау-Дюрсо. Новороссийск, Россия); международная молодёжная конференция-школа «Современные проблемы прикладной математики и информатики» (г. Дубна, 2012 г.); XIII всероссийская конференция молодых ученых по математическому моделированию и информационным технологиям (г. Новосибирск 2012 г.); XVII молодежная научная конференция Объединения молодых ученых и специалистов (ОМУС-2013) (г. Дубна 2013 г.); международная конференция по

математической теории управления и механике (г. Суздаль, Россия, 5-9 июля 2013 г.); The International Conference MATHEMATICAL MODELING AND COMPUTATIONAL PHYSICS (MMCP 2013) (г. Дубна 2013 г.).

Доклад на XVI международной конференции студентов, аспирантов и молодых учёных «Ломоносов» (г. Москва 2009 г.) отмечен грамотой за лучший доклад. Доклад на XVII молодежной научной конференции Объединения молодых ученых и специалистов (ОМУС-2013) (г. Дубна 2013 г.) отмечен почетным дипломом за лучший доклад секции «Математическое моделировании и вычислительная физика».

Публикации. Материалы диссертации опубликованы в 19 печатных работах, из них 4 статьи в журналах из перечня ведущих рецензируемых научных журналов и изданий, рекомендованных ВАК, [1-4], 4 публикации в других научных изданиях [5-8] и 11 тезисов докладов.

Личный вклад автора. Содержание диссертации и основные положения, выносимые на защиту, отражают персональный вклад автора в опубликованные работы. Подготовка к публикации полученных результатов проводилась совместно с соавторами в работах [1, 2, 4— 10]. Автору принадлежит основной вклад в разработку математической модели и алгоритмов выделения особенностей в полях газодинамических функций. Комплекс программ для выделения разрывов в различных газодинамических течениях был полностью реализован лично автором. Часть расчетов газодинамических течений, использованных в работе, также выполнена автором.

Структура и объем диссертации. Диссертация состоит из введения, четырех глав, заключения, библиографии и приложений. Общий объем диссертации 125 страницы, включая 66 рисунков и 1 таблицу. Библиография включает 48 наименований на 4 страницах.

Содержание работы

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

В первой главе приведен обзор известных методов применения вейвлет - анализа для выделения особенностей. Некоторые из этих методов (параграфы 1.3, 1.4) имеют достаточное

теоретическое обоснование, но не могут быть использованы на практике, так как либо требуют трудоемких расчетов (вычисление непрерывного вейвлет - преобразования), либо не приспособлены для анализа дискретных численных полей, с которым приходится иметь дело в большинстве прикладных задач. На практике могут быть использованы только некоторые дискретные модификации этих алгоритмов. Также приведено описание и сравнение некоторых прикладных алгоритмов (параграфы 1.5, 1.6), включая многомасштабную версию известного детектора Кэнни. В параграфе 1.7 описана модификация базиса вейвлетов, которая лучше адаптирована для представления двумерного сигнала и может быть использована для построения более эффективных, но и более ресурсоемких, алгоритмов выделения и анализа особенностей. В параграфе 1.8 приведен метод, предложенный С.Б. Базаровым, который не использует вейвлеты, этот метод был реализован в ходе исследований и апробирован в главе 2 (пункт 2.3.7). В параграфе 1.9 представлены методы локализации особенностей газодинамических полей, разработанные на первом этапе исследования данной задачи.

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

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

Если рассмотреть простейшие примеры разрывов (рисунок 1), сглаженных за счет свертки с гауссианом, а также три производных от сглаженных функций, видно, что сильному разрыву соответствуют нуль второй производной и локальный максимум модуля третьей производной. Слабому разрыву наоборот соответствуют нуль третьей производной и локальный максимум модуля второй производной. Для вычисления сглаженных второй и третьей производных от исходного поля могут быть использованы следующие выражения: '"„ =аг„ с„ • где ~~ исходное поле, а -

к т к т

постоянный коэффициент, {¿¡} и (г/,} - вещественная и мнимая части низкочастотного фильтра комплексного вейвлета Добешиб (цифра 6 обозначает, что фильтры вейвлета имеют 6 ненулевых элементов), {/•,.} и {} — сглаженные вторая и третья производные. {с, } -высокочастотный фильтр классического вещественного вейвлета Добешиб.

8

Сильный разрыв Слабый разрыв

Рисунок 1 — Шаблоны разрывов

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

Из соотношений Гюгонио на разрыве следует, что на контактных разрывах существует только разрыв плотности, но нет разрыва давления. Но поскольку существует некоторая погрешность при локализации разрывов, удобно определить новые операции над множествами точек, являющимися подмножествами узлов расчетной сетки. Введем натуральный параметр /?, характеризующий допустимую погрешность локализации. Расстоянием между узлами сетки и хл будем называть

...../.) = «Н* —.....К

Для операции пересечения (С = А П В ): точка аеС, если а е А и 36 е 8, такое что расстояние между а и Ь меньше К. Для операции разности (С = А \ В): точка аеС, если а е А и 3)е8, такого что расстояние между а и Ь меньше Л.

Рисунок 2 - Схема классификации узлов сетки

Ударным волнам соответствует множество точек: (\iZDCiCED) Т\{М2Р(~\СЕР). Так как на контактных разрывах рвется плотность, но не рвется давление, им соответствует множество точек: (Л/2£>ПС£0) \ (\tZPHCEP). Волнам разрежения соответствует множество точек: (МЕОГ\СгО)П(МЕРПСгР).

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

Результаты второй главы опубликованы в работах [1. 5].

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

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

ю

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

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

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

Установлено, что детектор чувствителен к дефектам алгоритмов расчета, благодаря чему может использоваться для оценки его корректности и качества.

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

Результаты третьей главы опубликованы в работах [3, 6].

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

Процедуру выделения и классификации разрывов в расчете можно разбить на четыре

этапа:

1) разделение расчетной сетки на ломаные,

2) обработка ломаных с помощью вейвлетов и выделение особенностей,

3) объединение результатов обработки ломаных,

4) фильтрация артефактов и классификация особенностей.

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

На втором этапе в соответствии с моделью выделения разрывов, предложенной во второй главе, производится независимая обработка ломаных. Обрабатываются только ломаные, у которых число узлов N больше 6, это число определяется количеством ненулевых элементов в фильтрах вейвлетов, используемых при обработке. Массивы плотности с1к и давления рк, где 0 < А: < V, заданные в узлах ломаной, также обрабатываются независимо. Для каждого из массивов вычисляется два преобразования т(х) и с(х)\

">М ' =Е<Г'-«Х#„-2Л,ДЛЯ 0 < к < V.

Н I и /

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

1) 'нули' - если значения массива в двух соседних узлах имеют разный знак или

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

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

Из множества нулей исключаются точки, соответствующие осцилляциям. Считается, что точка соответствует осцилляциям, если слева и справа от нее в радиусе трех точек есть выделенные нули. Таким образом, каждый узел ломаной получает некоторый набор из 8 возможных меток: mzd, ced, mzp, сер. med, ezd, шер, czp. Символ 'm' означает, что метка относится к основному детектору, 'с' — к корректору, 'd' означает, что метка характеризует поле плотности, а 'р' - давления. Символы 'z' и 'е' определяют, какие структуры были выделены в детекторе: 'z' соответствует переходам через ноль, а 'е' - локальным экстремумам модуля.

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

Процесс удаления артефактов основывается на предположении, что разрывы представляют собой набор кусочно-гладких кривых (поверхностей). Классификация разрывов производится на основании модели введенной во второй главе.

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

1) невозмущенный неподвижный газ: jc > 2, pi=l, Ui=0, pi=l,

2) объемная часть разряда: 0 < jc < 2, pj=l, чг=0, р2=12.4625,

3) область за фронтом падающей волны: х < 0, рэ=3.7629, Uj=2.5194, рз=9.6450.

Рисунок 3 - Распределение плотности в расчете на равномерной подвижной сетке, содержащей 2048 ячеек (слева), и разрывы, локализованные в расчете (справа).

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

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

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

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

Использование адаптивного подхода позволило существенно повысить качество расчета за счет уменьшения зон размазывания разрывов (особенно контактных рисунок 4).

Рисунок 4 - Распределение плотности в момент / =5 для адаптивного (черный) и неадаптивного (серый) расчетов.

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

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

Выполнено моделирование тесений для задач задача о сверхзвуковом обтекании осесимметричных тел. Были проведены расчеты, соответствующие уравнениям Эйлера и Рейнольдса с к — е моделью турбулентности. Параметры набегающего потока в обоих случаях одинаковы и соответствуют числу Маха 1.5. Число Рейнольдса для вязкой среды приближенно равно 3 ООО ООО. Расчеты выполнялись с помощью пакета №С2КЕ на треугольной сетке, в области с криволинейными границами.

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

структура разрывов за обтекаемым телом.

Рисунок 5 - Разрывы (сверху), локализованные в расчетах, выполненных по идеальной (слева) и вязкой (справа) моделям, и распределение градиента плотности (снизу)

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

При наложении результатов локализации структур в двух расчетах видно, что положение головной ударной волны в расчетах в точности совпадает. Также в обоих расчетах присутствует разрыв АВ, порожденный взаимодействием головной ударной воны с границей области и разрыв СО.

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

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

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

Результаты четвертой главы опубликованы в работах [2. 4, 7, В].

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

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

Публикации в журналах из перечня ведущих рецензируемых научных журналов и изданий, рекомендованных ВАК

1. Плёнкин, А. В. Локализация разрывов в полях газодинамических функций с помощью вейвлет анализа / А. Л. Афендиков, Л. И. Левкович-Маслюк, А. Е. Луцкий, А. В. Плёнкин // Математическое Моделирование. - 2008. - № 7. - С. 65-84.

2. Плёнкин, А. В. Вейвлетньш анализ локализованных структур в идеальной и вязкой моделях / А. Л. Афендиков, А. Е. Луцкий, А. В. Плёнкин // Математическое Моделирование. - 2011. -№ 1. - С. 41-50.

3. Плёнкин, А. В. Кратно-масштабный анализ газодинамических полей / А. В. Плёнкин // Вестник Московского университета Серия 1. Математика. Механика. -2011. -№ 2. -С. 56-59.

4. Плёнкин, А. В. Локализация особенностей газодинамических полей и адаптация расчетной сетки к положению разрывов / А. Л. Афендиков, А. Е. Луцкий, А. В. Плёнкин // Математическое Моделирование. - 2012. - № 12. - С. 49-54.

Публикации в других научных изданиях

5. Плёнкин, А. В. Локализация сингулярносгей газодинамических полей при помощи комплексных и вещественных вейвлетов / А. Л. Афендиков, В. В. Горбунова, Л. И. Левкович-Маслюк, А. В. Плёнкин // Препринт Института прикладной математики имени М. В. Келдыша РАН. - 2005. - № 98. - 32 с.

6. Плёнкин, А. В. Многомасштабный анализ особенностей газодинамических полей / А. Л. Афендиков, А. Е. Луцкий, А. В. Плёнкин // Препринт Института прикладной математики имени М. В. Келдыша РАН. - 2008. - № 98. - 17 с.

7. Плёнкин, А. В. Локализованные структуры в идеальной и вязкой моделях. Вейвлетньш анализ / А. Л. Афендиков, А. Е. Луцкий, А. В. Плёнкин // Препринт Института прикладной математики имени М. В. Келдыша РАН. -2009. -№ 78. - 10 с.

8. Плёнкин, А. В. Применение вейвлет анализа для выделения структур в расчетах газодинамических течений и для адаптации сеток / А. Л. Афендиков, А. Е. Луцкий, А. В. Плёнкин // Научная визуализация. - 2012. - № 3. - С. 8-25.

Напечатано в Типографии на Брестской ООО "Диалог", 123056, г. Москва, ул. 2-я Брестская, д. 39, стр. 3, ИНН: 7705928320 КПП: 771001001 Подписано в печать «17» марта 2014 г. Формата А5 (148x210мм)

Тираж 100 экз. Заказ №818 от 17.03.14

МОСКОВСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ ИМЕНИ М. В. ЛОМОНОСОВА

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

ПЛЁНКИН АНДРЕЙ ВАЛЕРЬЕВИЧ

04201457091

РАЗРЫВЫ ГАЗОДИНАМИЧЕСКИХ ФУНКЦИЙ В МЕТОДАХ СКВОЗНОГО СЧЕТА, ИХ АЛГОРИТМИЧЕСКАЯ ЛОКАЛИЗАЦИЯ И

КЛАССИФИКАЦИЯ

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

комплексы программ

ДИССЕРТАЦИЯ на соискание ученой степени кандидата физико-математических наук

Научный руководитель - профессор Афендиков Андрей Леонидович

Москва

2013

Оглавление

Введение..................................................................................................................................................4

1 Обзор известных методов выделения особенностей................................................................10

1.1 Основные определения и некоторые свойства вейвлетов................................................11

1.1.1 Показатель гладкости Липшица.................................................................................11

1.1.2 Условие Фурье..............................................................................................................12

1.1.3 Вейвлеты.......................................................................................................................13

1.1.4 Нулевые моменты вейвлетов......................................................................................14

1.1.5 Многомасштабный дифференциальный оператор....................................................15

1.1.6 Многомасштабный (кратномасштабный) анализ.....................................................16

1.1.7 Масштабирующая (скейлинг) функция и материнский вейвлет.............................17

1.1.8 Вейвлет - фильтры и дискретное вейвлет - преобразование...................................18

1.2 Связь асимптотики убывания вейвлет - преобразования с равномерной гладкостью Липшица на отрезке.........................................................................................................................19

1.3 Точечная гладкость Липшица.............................................................................................21

1.3.1 Теорема Жаффара........................................................................................................21

1.3.2 Конус влияния..............................................................................................................22

1.4 Максимумы модуля вейвлет - преобразования................................................................23

1.4.1 Максимумы модуля и выделение особенностей.......................................................23

1.4.2 Распространение максимумов.....................................................................................26

1.5 Многомасштабное выделение перепадов..........................................................................26

1.5.1 Алгоритм Кэнни выделения перепадов.....................................................................27

1.5.2 Многомасштабное выделение перепадов..................................................................28

1.5.3 Кривые максимумов.....................................................................................................29

1.5.4 Гладкость Липшица.....................................................................................................30

1.6 Увеличение четкости изображений, при помощи Лапласиана........................................32

1.7 Риджлеты...............................................................................................................................34

1.7.1 Оценки аппроксимации функции в различных базисах...........................................34

1.7.2 Базис риджлетов...........................................................................................................34

1.8 Выделение особенностей методом С.Б. Базарова.............................................................35

1.9 Локализация сингулярностей газодинамических полей при помощи вейвлетов..........36

1.9.1 Детектор, основанный на оценке показателя гладкости Липшица в точках исследуемого поля.......................................................................................................................37

1.9.2 Детектор на основе пары вещественных ортогональных вейвлет - фильтров.......39

1.10 Заключение...........................................................................................................................45

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

2.1 Введение................................................................................................................................46

2.2 Применение вейвлетов для выделения и классификации особенностей газодинамических полей.................................................................................................................47

2.2.1 Детектор на основе симметричных комплексных вейвлетов Добеши...................48

2.2.1.1 Симметричные комплексные вейвлеты Добеши..................................................48

2.2.1.2 Конструкция детектора............................................................................................49

2.2.1.3 Многомерная версия детектора..............................................................................52

2.2.2 Построение детектора-корректора.............................................................................54

2.2.3 Классификация разрывов............................................................................................55

2.2.3.1 Определения и обозначения....................................................................................55

2.2.3.2 Классификация.........................................................................................................56

2.2.4 Алгоритм локализации и классификации разрывов в двумерном случае, оценка числа необходимых вычислительных операций.......................................................................57

2.3 Численные эксперименты...................................................................................................59

2.3.1 Тестовый пример..........................................................................................................59

2.3.2 Численное моделирование..........................................................................................60

2.3.3 Исследование одномерной задачи..............................................................................61

2.3.4 Исследование 2Т> задачи..............................................................................................63

2.3.5 Выделение и классификация разрывов......................................................................64

2.3.6 Сравнение результатов локализации разрывов в расчетах, полученных по схемам 1-го и 2-го порядков аппроксимации.........................................................................................67

2.3.7 Обработка исходных данных методом С.Б. Базарова..............................................69

2.4 Заключение...........................................................................................................................71

3 Многомасштабный анализ особенностей газодинамических полей и оценка качества работы детектора..................................................................................................................................72

3.1 Введение................................................................................................................................72

3.2 Численное моделирование..................................................................................................73

3.3 Многомасштабный вейвлет анализ....................................................................................74

3.4 Структуры, локализуемые у границы расчетной области................................................83

3.5 Локализация разрывов с помощью корректора................................................................83

3.6 Сравнение расчетов выполненных на различных сетках.................................................86

3.7 Заключение...........................................................................................................................88

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

4.1 Введение................................................................................................................................90

4.2 Алгоритм...............................................................................................................................91

4.3 Численное моделирование..................................................................................................95

4.4 Анализ расчетов...................................................................................................................96

4.5 Локализация разрывов в трехмерных расчетах...............................................................101

4.6 Применение детектора для адаптации расчетной сетки к положению разрывов........105

4.7 Заключение.........................................................................................................................110

Заключение.........................................................................................................................................112

5 Приложения................................................................................................................................115

5.1 Симметричные комплексные вейвлеты Добеши............................................................115

5.2 Связь между вещественной и мнимой частями скейлинг функции..............................118

5.3 Лемма о сходимости масштабирующих функций..........................................................119

5.4 Низкочастотные фильтры использованных в работе вейвлетов...................................120

Список литературы............................................................................................................................122

Введение

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

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

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

Можно выделить две существенных особенности методов сквозного счета:

1) размазываются разрывы газодинамических функций;

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

Первая из указанных выше проблем может быть решена за счет измельчения расчетной сетки. Однако увеличение точности расчета только за счет равномерного измельчения сетки не оптимально [2, 3] и не всегда возможно даже на современных супер ЭВМ. В такой ситуации оказалось эффективным использование адаптивных сеток сгущающихся в окрестности разрывов. Эффективность применения адаптивных сеток при расчете газодинамических течений продемонстрирована в работах [4, 5, 6]. Использование адаптивных сеток, сгущающихся в областях высоких градиентов, позволило сократить число расчетных узлов в 30 раз, по сравнению с рав-

номерной сеткой, без ухудшения качества расчета.

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

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

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

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

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

К нашей постановке задачи ближе всего лежат исследования H.H. Яненко, Е.В. Ворож-цова[7] и С.Б. Базарова [8]. Отметим, что имеется значительное количество работ, посвященных решению задачи выделения краев (edge detection), возникающей при обработке изображений. Показательным является то, что в основе большинства этих методов лежит использование вейвлет - анализа, который в задачах обработки изображений доказал свою эффективность. Использование вейвлет - разложения также оказалось эффективным средством для ускорения расчетов использующих методы типа Галеркина [9].

Представленная работа состоит из четырех глав.

Глава 1 посвящена обзору известных методов выделения особенностей, в основе которых лежит применение вейвлет - анализа.

Глава 2 посвящена построению математической модели алгоритма выделения и классификации разрывов для одномерных течений газа. Модель основана на использовании построенного в [10] семейства комплексных симметричных ортогональных вейвлетов Добеши. Свойства этих вейвлетов (приложения 5.1, 5.2) оказались полезным инструментом в исследуемой задаче. Предобработка изучаемых полей использует соображения близкие к применяемым в [10] для повышения контрастности рентгеновских снимков.

Существенным отличием предлагаемого анализа от метода, предложенного в [10], является использование корректора, построенного на основе классических несимметричных вейвлетов и дополнительное использование соотношений на разрывах, наличие которых жестко связывает особенности различных компонент вектора газодинамических величин.

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

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

Альтернативный подход предложен в [11, 12, 13]. Там информация об особенностях течений извлекалась с помощью кратномасштабного анализа из коэффициентов вейвлет - разложения изучаемого поля по базису из комплексных ортогональных вейвлетов, а в качестве индикаторов особенностей использовались скачки фазы вейвлет - коэффициентов.

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

Отметим, что в расчетах могут существовать эффекты, соответствующие разным мае-

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

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

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

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