Инверсия данных электромагнитного каротажа на графических ускорителях

Вычислительные технологии

Том 20, № 1, 2015

Инверсия данных электромагнитного каротажа на графических ускорителях

В.Н. Глинских1'2'*, В. А. Горбатенко1

1 Новосибирский государственный университет, Россия

2Институт нефтегазовой геологии и геофизики им. А.А. Трофимука СО РАН, Новосибирск, Россия

*Контактный e-mail: GlinskikhVN@ipgg.sbras.ru

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

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

Введение

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

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

© ИВТ СО РАН, 2015

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

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

Для повышения быстродействия и производительности вычислений широко используются графические ускорители (ГПУ). Наиболее распространена технология NVIDIA CUDA (Compute Unified Device Architecture), предназначенная для проведения массивно-параллельных вычислений на ГПУ [4,5]. Основная идея этой технологии — представление вычислений в блочно-сеточном виде, что способствует гибкой адаптации алгоритмов на основе матрично-векторных вычислений. Технология даёт возможность использовать разные архитектурные особенности ГПУ, например, различные виды памяти и кэшей, эффективное применение которых позволяет существенно увеличить быстродействие и производительность вычислений. Численные решения задач скважинной геоэлектрики с использованием ГПУ рассматриваются в [6,7].

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

1. О подходах к решениям прямой и обратной задач электромагнитного каротажа

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

сти а (УЭП).

Каротажный прибор

/

Прискважинная зона

Пласт

Рис. 1. Интерпретационная геоэлектрическая модель

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

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

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

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

(•ехр ("Буи • г

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

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

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

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

5f = Абс,

Р({, g) < 1,

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

2. Особенности реализации и оптимизации вычислительного алгоритма

Технология NVIDIA CUDA предназначена для реализации вычислительных алгоритмов, предусматривающих параллелизм по данным. Алгоритм построения областей квазирешений обратной задачи с использованием полного перебора модельных параметров сведён к параллелизму по прямым задачам, решение которых выполняется независимо и одновременно каждым потоком. Для определения соответствия между исполнительными потоками используется уникальный идентификатор потока. Основная вычислительная сложность реализации решения прямой задачи — операция произведения матрицы чувствительностей на соответствующий вектор модельных параметров в соответствии с (1). Вычисление значения целевой функции (2) не является трудоёмким и предоставляется каждому потоку. В рамках рассматриваемой реализации алгоритма для хранения всех данных и результатов вычислений используется глобальная память. Далее рассмотрим несколько версий вычислительного алгоритма, реализованных при оптимизациях, направленных на сокрытие латентности глобальной памяти путём использования других типов памяти для различных видов данных.

Для анализа эффективности работы приложения на ГПУ применяется профилировщик платформы для разработчиков NVIDIA Nsight [15]. Здесь и далее результаты профилировки оптимизированных версий вычислительного алгоритма анализируются на модельном примере, в котором осуществляется перебор четырёх параметров при вариации каждого 50 раз и соответственно выполняется 6.25 млн решений прямых задач. Проводится сопоставление производительности вычислений, объёмов "чтения/записи" данных и пропускной способности при различных подходах к оптимизации алгоритма. Под объёмом "чтения/записи" данных понимается общий объём транзакций, а под скоростью — отношение данного объёма ко времени работы функции ядра на ГПУ. Приведённые ниже оценки этих показателей получены на оборудовании, включающем ЦПУ Intel Core Quad Q9300, оперативную память 4 Гб и графический ускоритель NVIDIA GTX Titan на чипе GK110.

Основная проблема при реализации вычислительного алгоритма — неоправданно большое количество операций "чтение/запись" с использованием глобальной памяти ГПУ. В результате профилировки алгоритма выяснено, что в рассматриваемом модельном примере производительность достигает 15 ГФлопс, а попадания L2 кэша, используемого только при обращениях мультипроцессора в глобальную память, составляют 44.8 %. Кроме регистровой другие виды памяти не применяются. При работе приложения происходит запись в глобальную память данных объёмом 82 Гб со скоростью 319 Гб/с, а чтение объёмом 32 Гб со скоростью 123 Гб/с.

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

При этом достигается производительность 90 ГФлопс, а попадания L2 кэша составляют 99.8 %. Оптимизация алгоритма приводит к существенным изменениям объёмов "чтения/записи", а также пропускных способностей. Происходит запись данных объёмом 440 Мб со скоростью 800 Мб/с, а чтение объёмом 264 Гб со скоростью 478 Гб/с.

С учётом особенностей разрабатываемого алгоритма и архитектуры ГПУ такие данные, как параметры перебираемых геоэлектрических моделей, после соответствующих вычислений каждым из потоков для хранения дублируются и в глобальную, и в разделяемую память. Так, для повышения быстродействия все обращения направляются к копии данных, находящихся в разделяемой памяти. Это решение обусловлено тем, что разделяемая память физически находится ближе к исполняющим элементам и доступ к ней относительно глобальной памяти осуществляется значительно быстрее. Поскольку в результате вычислений необходимо копировать все значения вектора модельных параметров, то копия, создаваемая в глобальной памяти, является также необходимой, поскольку разделяемая память доступна на "чтение/запись" только с ГПУ и копировать из неё в память ЦПУ невозможно. Благодаря такой оптимизации удалось достичь производительности 340 ГФлопс. Происходит запись в глобальную память данных того же объёма 440 Мб, а чтение уменьшается до 53 Гб. При этом из разделяемой памяти осуществляется чтение данных объёмом 141 Гб со скоростью 1 Тб/с.

При вычислениях используются некоторые исходные данные, а именно: диапазон и число вариаций модельных параметров, определяющие условия перебора. Согласованный доступ к указанным данным в глобальной памяти не может быть эффективно использован в силу существующего представления алгоритма, поэтому применяется текстурный кэш. Так как объём данных, хранимых в текстурной памяти, небольшой, то удаётся достичь 100%-го попадания текстурного кэша, благодаря чему скрывается время ожидания отклика памяти. В результате наблюдается увеличение производительности до 385 ГФлопс. В рассматриваемом модельном примере общий объём данных, считываемый из текстур, составляет 525 Мб.

Как уже отмечалось, вычислительный алгоритм с использованием полного перебора модельных параметров обладает высокой ресурсоёмкостью. Даже в простом модельном примере с перебором трёх параметров при реализации вычислений требуется копирование данных объёмом около 100 Мб из памяти ГПУ в память ЦПУ. Поэтому оптимизация этих транзакций является важной. Для увеличения скорости копирования применяется pinned-память на ЦПУ (page-locked память), особенность которой состоит в запрете на страничное хранение. С её использованием удаётся увеличить скорость копирования из памяти ЦПУ в память ГПУ с 1500 до 5000 Мб/с, а в обратном направлении — с 1700 до 5650 Мб/с. В силу больших объёмов данных время их копирования становится значительным. Для уменьшения этих временных затрат используется подход, заключающийся в одновременном вычислении и копировании уже полученных данных. Для организации вычислений и копирования с перекрытием необходимо общий объём вычислений разделить на независимые части и вызывать их последовательно в разных потоках CUDA (CUDA Stream), используя при этом асинхронное копирование в том же потоке CUDA (рис. 2). Приведённый рисунок построен с применением инструмента профилировки NVIDIA Visual Profiler [16].

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

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

В [0] GeForce GTXTTTAN В Context 1

МегпСру (HtoD) L V МегпСру (DtoH) Щ Compute В Streams L Default

calcF kerfrnt, rnt, rnt, rnt, rnt, rnt*, float* float*, float*, flc...

calcF.kerCint, rr>t_ int, ir>tH rrrtH int* float* float*, float* flo..

B [0] GeForce GTXTTTAN B Context 1

B Compute

B Streams L Default L Stream 13 L Stream L4 L Stream IS L Stream 16

calch_kertin... calcF_ker(int...

carcF_ker(rnt... calcF_ker[rnt..

calcF_ker(in,,,

calcF_ker(int...

cafcF_ker(int...

carcF_ker(int...

Рис. 2. Организация вычислений и копирования данных без перекрытия (а) и с перекрытием (б) на ГПУ

2000

о

1500

к К

<D

Ч

О К

tr

PQ

1930

1000

CD Он

PQ

500

363

1 114 112 99 70

0 1 2 3 4 5

Номер оптимизации вычислительного алгоритма

а

Рис. 3. Быстродействие версий вычислительного алгоритма при оптимизации на ГПУ

удаётся достичь скорости 950 Гб/с, что почти в 3 раза превышает этот показатель при использовании глобальной памяти.

Оценки быстродействия описанных выше оптимизаций вычислительного параллельного алгоритма на ГПУ представлены на рис. 3. Здесь время вычислений в столбце с номером 0 соответствует версии алгоритма без оптимизации, 1 — версии, в которой происходит отказ от хранения синтетических данных в глобальной памяти, столбец 2 — в качестве кэша используется разделяемая память, 3 — текстурный кэш применяется для хранения небольшого объёма данных, 4 — проводится оптимизация копирования данных с перекрытием вычислений, 5 — текстурный кэш вовлекается для обращений к матрице чувствительностей.

3. Оценки быстродействия и производительности вычислений

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

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

Для сравнительного анализа при вычислениях на ЦПУ используется аналогичный описанному выше последовательный алгоритм. Этот алгоритм не оптимизирован и применяется для некоторой общей оценки "сверху" эффективности оптимизированных версий реализованного параллельного алгоритма при вычислениях на ГПУ.

В таблице приведены оценки быстродействия и производительности вычислительного алгоритма в зависимости от числа варьируемых модельных параметров. На рис. 4 представлены основные характеристики эффективности вычислительного алгоритма численной инверсии при определении трёх, шести и девяти модельных параметров. Здесь показаны времена (а), производительности и ускорения (б) вычислений на ГПУ и ЦПУ в зависимости от числа варьируемых модельных параметров. Из приведённых в таблице и на рисунке данных следует, что время вычислений на ГПУ существенно ниже чем, на ЦПУ. Ускорение возрастает с увеличением числа модельных параметров и достигает наилучших показателей при вариациях семи параметров. Существенный рост производительности и ускорения связан с тем, что современные ГПУ содержат тысячи потоковых процессоров и для полной загрузки устройства требуются сотни тысяч потоков. Последующее снижение производительности и соответственно ускорения определяется увеличением вычислений и повышением размерности задачи.

Быстродействие и производительность вычислительного алгоритма в зависимости от числа

варьируемых модельных параметров

Характеристика вычислений Число варьируемых модельных параметров

3 4 5 6 7 8 9

Время вычислений на ГПУ, мс 1.5 1.6 5.7 6.2  101 5.1  102 6.1  103 7.1  104

Время вычислений на ЦПУ, мс 4.6  101 5.3  102 5.8  103 6.7  104 7.2  105 7.9  106 8.5  107

Производительность вычислений, ГФлопс 5.7 53.0 277.2 364.2 513.7 421.4 390.2

Быстродействие (ускорение), раз 30 354 1017 1079 1408 1293 1204

«

s д

<L>

ч о

s

F PQ

a PQ

108 106 104 102 10°

ГПУ  ЦПУ

1 1

1 il _

3 6 9

Число варьируемых параметров

1200

900

«

s к

CD

4 о

5 F

PQ

§ 600 s

H о

5

cd

6

X

1 1  Производительность, ГФлопс  Ускорение, раз Ц

I

il _

3 6 9

Число варьируемых параметров

Рис. 4. Основные характеристики эффективности вычислительных алгоритмов: а — время вычислений параллельного алгоритма на ГПУ и последовательного на ЦПУ, б — производительность и ускорение параллельного относительно последовательного

400

>g 300

s

х

<D

ч

0

Е 200

hP

PQ g

1 100

PQ

GT200

_ _ GF110 GK104

_ GK110

280 560 680 Titan Графические ускорители NVIDIA GTX

Рис. 5. Быстродействие вычислительного алгоритма на ГПУ с различной архитектурой

а

Далее был выполнен сравнительный анализ быстродействия реализованного параллельного алгоритма на ГПУ в модели флюидонасыщенного коллектора при определении девяти параметров с построением областей квазирешений обратной задачи на ГПУ разных архитектур и поколений (рис. 5). Вычислительные тесты проведены как на современных графических ускорителях NVIDIA GTX Titan на основе процессора GK110 и NVIDIA GTX680 с чипом GK104, имеющих архитектуру Kepler, так и на графических ускорителях предыдущих поколений NVIDIA GTX560 (чип GF104, архитектура Fermi) и NVIDIA GTX280 (чип GT200, архитектура Tesla). Как видно из рисунка, увеличение быстродействия вычислений на GTX Titan по сравнению с GTX280 составляет около 5 раз. Это обусловлено интенсивным развитием графических ускорителей, главным образом заключающимся в значительном наращивании числа вычислительных ядер с понижением их тактовой частоты.

4. Результаты численной инверсии на интервале нефтяных пластов

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

По зашумленным синтетическим данным электромагнитного каротажа с учётом (3) построены и проанализированы области квазирешений обратной задачи в модели флю-идонасыщенного коллектора с переменным по глубине нефтесодержанием и сложным строением прискважинной зоны (рис. 6). Для этого примера приведены интересные с практической точки зрения области квазирешений обратной задачи по УЭП зон про-

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

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

Заключение

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

Выявлены различные подходы к оптимизации параллельного алгоритма на графических ускорителях, основанные на детальном анализе производительности вычислений и эффективности работы с памятью, результатом которых стало сокрытие латентности глобальной памяти посредством использования разделяемой и текстурной. Получены оценки быстродействия и производительности вычислений оптимизированных версий параллельного алгоритма на графических ускорителях NVIDIA с архитектурами различных поколений (Kepler, Fermi, Tesla). Показано, что по сравнению с идентичным последовательным алгоритмом применение предложенного параллельного алгоритма на современных графических ускорителях позволяет достигать ускорения вычислений до 1400 раз, а производительности — до 530 ГФлопс.

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

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

[1] Табаровский Л.А., Эпов М.И., Сосунов О.Г. Оценка разрешающей способности электромагнитных методов и подавление помех в системах многократного наблюдения (теория, алгоритмы, программы) // Новосибирск, 1985 (Препр. АН СССР. Сиб. отд-ние. ИГиГ; № 7). 48 с.

Tabarovsky, L.A., Epov, M.I., Sosunov, O.G. Estimating resolving power of electromagnetic methods and noise reduction in systems with multiple observations (theory, algorithms and programs) // Novosibirsk: Preprint IGIG SO AN SSSR, 1985, No. 7. 48 p. (in Russian).

[2] Табаровский Л.А., Эпов М.И. Оценка разрешающей способности электромагнитных методов // Геология и геофизика. 2006. Т. 47, № 5. С. 568-578.

Tabarovsky, L.A., Epov, M.I. Estimating resolution of electromagnetic methods // Geology and Geophysics. 2006. Vol. 47, No. 5. P. 568-578. (in Russian).

[3] Эпов М.И., Глинских В.Н., Никитенко М.Н., Сухорукова К.В. Быстрое решение двумерной обратной задачи высокочастотного электромагнитного каротажа в нефтегазовых скважинах // Геофиз. журнал. 2012. Т. 34, № 4. С. 292-297.

Epov, M.I., Glinskikh, V.N., Nikitenko, M.N., Suhorukova, C.V. Fast solution of the two-dimensional high-frequency electromagnetic logging inverse problem in oil and gas wells // Geophysical Journal. 2012. Vol. 34, No. 4. P. 292-297 (in Russian).

[4] CUDA C Programming Guide. Design Guide (PG-02829-001_v6.0, February 2014) // NVIDIA CUDA Toolkit Documentation. 2014. Available at: http://docs.nvidia.com/cuda/ pdf/CUDA_C_Programming_Guide.pdf (date of visit: 01.09.14).

[5] Боресков А.В., Харламов А.А. Параллельные вычисления на GPU. Архитектура и программная модель CUDA. М.: Изд-во Московского ун-та, 2012. 332 с.

Boreskov, A.V., Kharlamov, A.A. Parallel Computations on GPU. CUDA Architecture and Software Model. Moscow: MSU, 2012. 332 p (in Russian).

[6] Глинских В.Н., Эпов М.И., Лабутин И.Б. Моделирование диаграмм электромагнитного каротажа на графических процессорах // Вычисл. технологии. 2008. Т. 13, № 6. С. 50-60.

Glinskikh, V.N., Epov, M.I., Labutin, I.B. Simulation of electromagnetic logs using graphics processing unit // Comput. Technologies. 2008. Vol. 13, No. 6. P. 50-60 (in Russian).

[7] Лабутин И.Б., Суродина И.В. Разработка параллельных алгоритмов для решения задач каротажа на графических процессорах // Параллельные вычислительные технологии (ПаВТ'2012): Труды международной научной конф. Челябинск: Издательский центр ЮУрГУ, 2012. С. 690-697.

Labutin, I.B., Surodina, I.V. Development of parallel algorithms for solving logging problems on GPUs // Parallel computational technologies 2012: transactions of the international scientific conference. Chelyabinsk: YuUrGU, 2012. P. 690-697 (in Russian).

[8] Суродина И.В., Эпов М.И. Влияние биополимерных буровых растворов на диаграммы высокочастотного электромагнитного каротажа // Геология и геофизика. 2012. Т. 53, № 8. С. 1062-1069.

Surodina, I.V., Epov, M.I. High-frequency induction data affected by biopolymer-based drilling fluids // Russian Geology and Geophysics. 2012. Vol. 53, No. 8. P. 817-822.

[9] Нечаев О.В., Шурина Э.П., Федорук М.П. Использование векторного метода конечных элементов для численного решения квазистационарных уравнений Максвелла // Вычисл. технологии. 2004. Т. 9, № 5. С. 73-81.

Nechaev, O.V., Shurina, E.P., Fedoruk, M.P. Using the vector finite element method for numerical solutions of the quasi-steady Maxwell's equations // Comput. Technologies. 2004. Vol. 9, No. 5. P. 73-81 (in Russian).

[10] Еремин В.Н., Нечаев О.В., Хаберхауэр ШШ. и др. Параллельная реализация математического моделирования процессов электромагнитного каротажа зондовым комплексом ВИКИЗ // Там же. 2007. Т. 12, № 6. С. 18-33.

Eremin, V.N., Nechaev, O.V., Haberhauer, S. Parallel realization of mathematical modelling of electromagnetic logging processes using VIKIZ probe complex // Comput. Technologies. 2007. Vol. 12, No. 6. P. 18-33 (in Russian).

[11] Шурина Э.П., ШШтабель Н.В. Анализ векторных конечноэлементных аппроксимаций уравнений Максвелла в анизотропных средах // Там же. 2013. Т. 18, № 4. С. 92-104. Shurina, E.P., Shtabel, N.V. Analysis of vector finite element approximations of Maxwell's equations in anisotropic media // Comput. Technologies. 2013. Vol. 18, No. 4. P. 92-104 (in Russian).

[12] Эпов М.И., Глинских В.Н. Электромагнитный каротаж: Моделирование и инверсия. Новосибирск: ГЕО, 2005. 98 с.

Epov, M.I., Glinskikh, V.N. Electromagnetic Logging: Simulation and Inversion. Novosibirsk: GEO, 2005. 98 p. (in Russian).

[13] Эпов М.И., Никитенко М.Н. Система одномерной интерпретации данных высокочастотных индукционных каротажных зондирований // Геология и геофизика. 1993. Т. 34, № 2. С. 124-130.

Epov, M.I., Nikitenko, M.N. System for one-dimensional interpretation of high-frequency induction logging sounding data // Geology and Geophysics. 1993. Vol. 34, No. 2. P. 124-130 (in Russian).

[14] Глинских В.Н., Никитенко М.Н., Эпов М.И. Линеаризованные решения двумерных прямой и обратной задач высокочастотного электромагнитного каротажа в проводящих средах с учётом токов смещения // Там же. 2013. Т. 54, № 12. С. 1942-1951. Glinskikh, V.N., Nikitenko, M.N., Epov, M.I. Processing high-frequency electromagnetic logs from conducting formations: linearized 2D forward and inverse solutions with regard to eddy currents // Russian Geology and Geophysics. 2013. Vol. 54, No. 12. P. 1515-1521.

[15] NVIDIA Nsight Visual Studio Edition 4.0 User Guide // NVIDIA CUDA Technical Documentation. 2014. Available at: http://docs.nvidia.com/gameworks/index.html#developer-tools/desktop/ nsight_visual_studio_edition_user_guide.htm (date of visit: 01.09.14).

[16] Profiler User's Guide (DU-05982-001_v6.0, February 2014) // NVIDIA CUDA Technical Documentation. 2014. Available at: http://docs.nvidia.com/cuda/pdf/CUDA_Profiler_Users_Gui-de.pdf (date of visit: 01.09.14).

Поступила в 'редакцию 4 сентября 2014 г.

Electromagnetic logging data inversion on GPU

Glinskikh, Viacheslav N.1,2'*, Gorbatenko, Vadim A.1

1 Novosibirsk State University, Novosibirsk, 630073, Russia 2Trofimuk Institute of Petroleum Geology and Geophysics SB RAS, Novosibirsk, 630090, Russia * Corresponding author: Glinskikh, Viacheslav N., e-mail: GlinskikhVN@ipgg.sbras.ru

The paper deals with the development of algorithms and software for processing of electrical logging data and interpretation on GPUs.

The traditional scheme for solving the inverse problem of borehole geoelectrics is based on the targeted search of model parameters. It includes experimental and simulated data matching, and utilizes iterative solution of the direct problem. Ambiguous relationship between electromagnetic field characteristics and geological medium properties along with inaccuracy of measured data leads to imprecise determination of parameters in the resulting model. As a result, instead of an unique solution there are a great number of equivalent models that can be described by quasi-solution areas in the inverse problem.

For creating such quasi-solution areas, we use a complete enumeration of model parameters. However, a sequential computing algorithm is utterly resource-intensive

© ICT SB RAS, 2015

but inefficient. Substantial improvement of computational efficiency may be achieved with the help of a parallel algorithm on multiprocessors.

We have implemented a parallel algorithm for the numerical solution of the electromagnetic logging of two-dimensional inverse problem on GPUs, which utilizes NVIDIA CUDA technology. Along with that, we have revealed various approaches to the optimization of the GPU parallel algorithm, based on detailed analysis of computational performance and efficiency for different memory types. The processing speed and performance of the optimized version of the parallel algorithm on NVIDIA GPUs with architectures of different generations (Kepler, Fermi, Tesla) are estimated. It is shown that implementation of this up-to-date GPU parallel algorithm speeds up calculations by a factor of 1400 and increases performance up to 530 Gflop/s in comparison with the identical sequential algorithm.

In addition, we have performed the numerical inversion of noisy synthetic data and carried out the comparative analysis of the inverse problem quasi-solution areas in complex-structure geologic models. Using the developed algorithm the electrical conductivities of invasion zones and non-invaded beds are reconstructed, which includes the estimation of inversion errors.

Keywords: parallel computations, GPU, NVIDIA CUDA technology, inverse problem, electromagnetic logging, electrical conductivity.

Received 4 September 2014