Курсовая работа (т): Моделирование распространения воздушных потоков и дыма с помощью графического процессора

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

Введение

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

Постепенно ситуация начала меняться и появились первые реализации фильтров пост-обработки, использующие GPU. Среди программного обеспечения, связанного с видео, такие фильтры одними из первых стали использовать графические ускорители. В силу естественности параллелизации GPU-реализации таких фильтров появились во всех профессиональных видеоредакторах, таких как, Adobe After Effects, Pinnacle Studio, Vegas.

Немного медленнее ситуация меняется среди программных систем, осуществляющих рендеринг, визуализацию готовых сцен. Первый коммерческий визуализатор, использующий GPU - NVIDIAOptiX - появился только в 2009 году, через 2 года после появления технологии CUDA, позволяющей проводить вычисления на GPU. Тем не менее на момент написания данного текста уже представлено несколько продуктов, способных осуществлять рендеринг, как полностью на GPU, так и в гибридном режиме - Chaos Group V-Ray RT, Octane, Arion, FurryballRT, Redshift, Cycles, moskito Render, Indigo Render, Thea Presto Render и т.д.

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

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

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

·        Чаще всего GPU используются только для отображении примерного вида спецэффекта в окне предпросмотра, а не для финального рендеринга.

·        Внешний вид эффектов, полученных при использовании GPU, отличается от вида эффектов, полученных с использованием CPU. В частности, финальная картинка, полученная на рендер-ферме с CPU слабо соответствует визуализации полученной в окне предпросмотра при использовании GPU.

При всем этом, несмотря на очевидные проблемы, художники по спецэффектам уже сейчас активно используют существующие плагины, например, ILM Plume, SideFX Houdini, Blender Cycles, поддерживающие видеоускорители, чтобы сэкономить время создания новой компьютерной Рисуноки. Еще большей пользы можно было бы добиться, если бы в широком доступе существовало ПО, способное производить не только предварительные расчеты, но и работать с большими сетками на GPU-кластерах, генерируя спецэффекты для финальной картинки. В моменты цейтнота и приближающегося срока сдачи того или иного фильма или ролика даже в крупных по размерах студиях художники подключают по сети свои персональные компьютеры к рендер-ферме, таким образом увеличивая суммарную производительность рендер-фермы и уменьшая временные затраты. При этом большое количество ресурсов в виде 3D-ускорителей простаивает вхолостую. Схожая ситуация наблюдается в маленьких ив домашних студиях. Там художник, обладая одним или несколькими мощными компьютерами лишь частично использует их ресурсы.

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

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

1. Обзор существующих работ

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

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

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

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

. Постановка задачи

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

где  и p - искомые вектор скорости и давление, а  - известная внешняя сила, например, описанная далее сила плавучести.

Так как плотность считается константной и, в данном случае, равной единице, то в рассматриваемой системе имеется четыре неизвестных и четыре уравнения. На внешних границах имеющей форму параллепипеда области Gext, в которой будут проводиться расчёты, заданы естественные условия, позволяющие газу беспрепятственно входить и выходить из области. На границах объектов произвольной формы, расположенных внутри области Gint наложены условия не протекания, препятствующие попаданию дыма внутрь объектов. Начальное распределение искомых величин задаётся аниматором исходя из того, что за эффект он хочет воспроизвести, и в общем случае представляется заранее известными функциями U0 и p0, которые должны быть согласованы с граничными условиями.

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

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

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

. Метод решения

Расщепление.

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

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

А в предварительный вектор скорости вносится соответствующая поправка:

В-третьих, необходимо найти новые значения для плотности и температуры дыма:

Адвекция.

Отдельно стоит остановиться на схеме численного решения уравнений адвекции, которые нужны для вычисления величин  Если для этого использовать аппроксимацию конечной разностью, как это сделано в выкладках выше, то не удастся обеспечить абсолютную устойчивость метода, а наличие данного свойства является крайне важным в силу ограниченности вычислительных мощностей. Поэтому в рассматриваемом алгоритме используется полулагранжевый подход, в рамках которого на этапе решения уравнений адвекции - и только на нём - жидкость представляется набором частиц, расположенных в центрах ячеек специальной сетки (MACgrid, см Рисунок 1) и обладающих соответствующими характеристиками.

Например, в случае температуры это означает, что для узла с номером (I, j, l) требуется вычислить величину:

,

где  - шаг по пространству. После этого значение на следующем временном слое будет определено как . Очевидно, что компоненты вектора  практически всегда будут дробными, поэтому для вычисления искомой величины температуры необходимо провести интерполяцию значений из восьми ячеек сетки, расстояние которых до узла (i −vx , j −vy , l -vz) не превышает. Если же индексы выходят запределы расчётной области, то соответствующее значение определяется исходя из граничных условий.

Рисунок 1. Одна ячейка MAC-сетки

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

Проекция.

Другим этапом решения исходной задачи является этап проекции. Он является наиболее ресурсоемким во всей предложенной схеме. Этот этап численно эквивалентен решению внешней задачи Неймана для уравнения Пуассона:

где:

Разностная аппроксимация данной задачи выглядит так:

А ее решение сводится к решению системы линейных алгебраических уравнений.

В данной работе не будет уделено внимание их реализации. Алгоритм работы этапа адвекции приведен лишь для введения понятия MAC-сетки, и описания требований к солверу для проекции. Основное внимание в этой работе будет уделено этапу проекции.

. Выбор метода выполнения этапа проекции

Для того чтобы удостовериться в осмысленности использования графических ускорителей, и чтобы выбрать наилучший численный метод решения системы, было решено сравнить методы решения СЛАУ из нескольких готовых библиотек. Для GPU была выбрана библиотека AmgX, а для CPU - Hypre. Эти библиотеки обладают наиболее полным набором возможностей для соответствующих платформа. В частности, каждая из них имеет реализацию алгебраического многосеточного метода (AMG) и стабилизированного метода бисопряженных градиентов (BiCGStab) - наиболее производительных устойчивых современных вычислительных методов.

Тестирование было произведено на 3 типах сцен: без объектов внутри, с расположенной по центру сферой, со 125 равномерно распределенными по сцене сферами. Для валидации использовалась следующая задача математической физики:

Рисунок 2. Время достижения требуемой невязки: AMG, пустая область.

Рисунок 3. Время достижения требуемой невязки: AmgX область с выколотой сферой

На машине с Intel Xeon E5620, 32 GB RAM и NvidiaTesla K40c все тесты показали превосходство AmgX над Hypre. Ускорение на GPU достигало 20 раз. Кроме того, было выявлено превосходство AMG над BiCGStab на всех невязках для всех типов задач.

Полный многосеточный метод.

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

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