Московский Государственный Университет им. М.В. Ломоносова

Факультет Вычислительной Математики и Кибернетики

Кафедра Суперкомпьютеров и Квантовой Информатики

Магистерская программа «Технологии параллельного

программирования и высокопроизводительные вычисления»

Магистерская диссертация

Исследование эффективности параллельных

алгоритмов моделирования распространения

акустических волн в гетерогенных средах в задачах

сейсморазведки для высокопроизводительных систем Работу выполнил студент Бурцев А. П.

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

к.ф-м.н. доц. Попова Н. Н., Ген. директор ООО «GEOLAB» Курин Е. А.

Москва Содержание Аннотация

Введение

Глава 1. Использование высокопроизводительных систем для моделирования геологических сред и анализа результатов сейсморазведки

1.1 Физические основы сейсморазведки

1.2 Обоснованность выбора технических средств

1.3 Обзор работ посвящённых моделированию геологических сред

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

1.5 Краткие выводы по главе

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

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

2.2 Структура параллельного алгоритма

2.3 Формат входных/выходных данных

2.4 Особенности программной реализации алгоритма на суперкомпьютере «Ломоносов»

2.5 Особенности программной реализации на GPU

2.6 Краткие выводы по главе

Глава 3. Исследование эффективности

3.1 Вычислительный эксперимент на базе IBM Blue Gene/P

3.2 Вычислительный эксперимент на базе ноутбуков ASUS N61VN и ASUS K53S.

3.3 Вычислительный эксперимент на базе суперкомпьютера «Ломоносов»........... 3.4 Краткие выводы по главе

Заключение

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

Дальнейшее развитие работы

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

Приложение А

Приложение Б

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

Проводится исследование эффективности полученных программ для большого числа вычислительных узлов. Осуществляется покомпонентная оценка элементов программ. Анализируются основные проблемы, возникающие при работе с большим числом вычислительных узлов, и предлагаются методы их решения. Программы тестировались на базе суперкомпьютеров «Ломоносов» и Blue Gene/P, а так же на графических ускорителях компании NVIDIA GT 240M и GT 610M.

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

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

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

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

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

OpenMP для SMP архитектуры, MPI для кластеров и CUDA для графических ускорителей компании NVIDIA [1].

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

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

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

В сейсморазведке различают два основные метода: метод отраженных волн (МОВ) и метод преломленных волн (МПВ). Меньшее применение находят методы, высокоточным определением геометрии геологического разреза (ошибки менее 1 %), стало возможным благодаря применению трудоемких систем возбуждения и наблюдения, обеспечивающих одновременный, иногда многократный съем информации с больших площадей и ее цифровую обработку на ЭВМ. Это обеспечивает выделение полезных, чаще однократно отраженных или преломленных волн среди множества волн-помех. Подобные задачи отличает большой объём обрабатываемых данных, что требует использовать не простые «домашние»

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

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

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

Рис. 1: Схема сейсморазведки на море. [3] Теория распространения упругих (сейсмических) волн базируется на теории упругости, так как геологические среды в первом приближении можно считать упругими. Поэтому напомним основные определения и законы теории упругости применительно к однородным изотропным средам.

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

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

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

При деформации в упругом теле возникают внутренние напряжения, обусловленные упругим взаимодействием между частицами тела. На каждую площадку малого размера, мысленно выделяемую в теле, действуют напряжения, имеющие в общем случае составляющую, перпендикулярную к площадке,— нормальное напряжение, и две, направленные вдоль площадки, называемые сдвиговыми напряжениями. Три компоненты напряжения задаются с помощью шести компонент тензора напряжения. Эти шесть компонент связаны с шестью компонентами малых деформаций законом Гука [4].

1.2 Обоснованность выбора технических средств.

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

Стоимость одной разведочной скважины на суше — 1 млн. долларов за Стоимость одной разведочной скважины на море — 10 млн. долларов за Средняя глубина скважины — 1.5 — 3 километра в настоящее время, но из-за истощения поверхностных запасов в ближайшем будущем эта величина составит 4 — 6 километров.

Высокие технологические риски, без предварительного изучения геологических данных велик риск разрушения не законченной скважины, а значит придётся начинать всё бурение сначала (См. Рис. 2).

Средняя стоимость высокопроизводительного вычислителя — 4 млн.

долларов (20Тфлопс, УГТУ, 2007г) Средняя стоимость часа работы — около 157 долларов (20Тфлопс, УГТУ, Рис. 2: Схема геологического участка на море. [3] Таким образом существует достаточное экономическое обоснование для использования высокопроизводительных вычислителей. Но не будет ли достаточно простых «персональных» машин?

Требования накладываемые реальной задачей (Карачаганак, См. Рис. 3) [3]:

Площадь наблюдений ~ 900 кв. км Количество пунктов возбуждения (ПВ) ~ Количество пунктов приёма (ПП) ~ Количество средних точек (ОСТ) ~ 9 Кратность сейсмограммы ОСТ Количество трасс (записей) ~ 3 миллиарда Длина записи 20 сек, 10000 отсчётов, 4 байта/отсчёт Объём данных > 100 терабайт Рис. 3: Месторождение Карачаганак, 3D съёмка, 2009г. [3] Решение задачи подобной сложности невозможно без использования специальных машин обладающих повышенной производительностью. Решение этого класса задач на «домашних» машинах невозможно.

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

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

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

Гибридный алгоритм предложенный в работе [6] также как и предыдущий использует специальную декомпозицию данных для достижения высокой эффективности в задачах имеющих не более 50% вязкоупругих включений. Кроме декомпозиции используются различные методы моделирования волн отдельно в вязкоупругих включениях и во всей остальной области. Алгоритм использует сравнительно небольшое число процессоров и основан на технологиях MPI и OpenMP.

Примером работы с уже существующими библиотеками (Intel®MKL) служит работа [7]. Важной особенностью предложенного подхода является использование разложение Лагерра по времени, вместо преобразования Фурье. Для организации параллельного счёта данные разбивается на прямоугольные области по которым счёт ведётся независимо, но возникает необходимость синхронизации пограничных зон.

Алгоритм использует возможности MPI для синхронизации и OpenMP для повышения скорости расчётов.

В работе [8] используется принцип геометрического параллелизма при распределении данных, а так же модель MPI+OpenMP, при моделировании распространения упругих волн в гетерогенных средах.

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

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

Примером полномасштабного исследования эффективности параллельного алгоритма сейсмической миграции в обратном времени для Blue Gene/P может служить работа [10]. В работе рассмотрены основные проблемы возникающие при использовании программы на разном числе процессоров, например такие как влияние кэша. Так же в работе проведено покомпонентное исследование элементов программы.

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

Задачи данной работы:

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

с учётом особенностей массивно-параллельной вычислительной системы с учётом особенностей суперкомпьютера «Ломоносов»

с учётом особенностей графических ускорителей распределённой памятью (Blue Gene/P и «Ломоносов»), при использовании достаточно большого числа вычислительных узлов, для определения ключевых проблем при организации масштабных вычислений.

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

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

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

x (горизонтальное удаление) и z (глубина), и шагом по времени t. Для граничного слоя реализованы отражающие условия — элементы поля давления всегда не возмущены (давление равно 0).

распространение акустической волны в твёрдом неоднородной среде [11], для двумерной нестационарной задачи с размерностью ( x, z, t) (1):

( x, z ) — плотность в каждой точке пространства.

v ( x, z ) — скорость акустической волны в каждой точке пространства.

момент времени t.

S ( x, z, t) — внешняя сила (функция источника) в точке пространства с координатами (x, z ), в момент времени t.

Для решения этой задачи воспользуемся методом конечных разностей заменив частные производные на конечно-разностные соотношения (Шаблон представлен на Рис. 4). Исходное уравнение примет вид (2):

Рис. 4: Графическая интерпретация конечно-разностной схемы используемой в данной работе.

Преобразовав последнее уравнение в явную конечно-разностную схему 2-ого порядка по времени получим (3):

излучения источников.

Заменив оператор Лапласа конечно-разностной схемой 4-ого порядка по пространственным координатам получим (4):

Скалярное произведение градиентов заменим конечной разностью (5):

Таким образом мы можем решать данную задачу имея информацию о двух начальных состояниях по времени и четырёх соседних точках по каждой из пространственных координат в момент времени состояния по времени:

граничные условия — во всех граничных точках поставленную задачу на ЭВМ.

Условием сходимости данной конечно-разностной схемы является критерий максимальный шаг по пространственным координатам: v =max ( v (x, z)) — максимальная скорость распространения акустической волны в среде.

последовательный вариант алгоритма, реализованного в рамках Open source проекта Madagascar [12].

Алгоритм состоит из следующих шагов:

1. Обработка и контроль корректности входных параметров.

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

3. Расчёт поля давлений (для каждого шага по времени):

1. Копирование данных, полученных на предыдущем шаге по времени, на 2. Расчёт новых значений поля давлений (текущий шаг по времени).

3. Обмен «теневыми» границами между соседними вычислительными В алгоритме используется геометрическая декомпозиция данных: каждый процессор обслуживает прямоугольную область. Нагрузка распределяется равномерно между всеми процессорами. Данные соседних процессоров пересекаются (См. Рис. 5).

Рис. 5: Схема распределения данных для четырёх процессоров.

Оценку оптимального числа вычислительных узлов следующим образом:

sizeof(float)( 3 m+ t+ Buffer sizem) Total available memory матрице, — число вычислительных узлов, оценку их максимального числа можно получить в результате пробного эксперимента на конкретной вычислительной вычислительной системе, из условия того, что время расчёта не должно превышать время обмена данными.

Все данные представляются в виде бинарных файлов. Файлы имеют следующий формат: число измерений (int), число элементов по каждой размерности (массив long int), число матриц (long int), матрицы представленные в виде линейного массива (массив float).