ё Д.А.Первухин, Ю.В.Ильюшин
Параллельный анализ геоданных гидролитосферных пластов.
УДК 681.5
ПАРАЛЛЕЛЬНЫЙ АНАЛИЗ ГЕОДАННЫХ ГИДРОЛИТОСФЕРНЫХ ПЛАСТОВ МИНЕРАЛЬНОЙ ВОДЫ КИСЛОВОДСКОГО МЕСТОРОЖДЕНИЯ НАРЗАНА
Д.А.ПЕРВУХИН, Ю. В. ИЛЬЮШИН
Санкт-Петербургский горный университет, Россия
Эколого-курортный район Кавказских Минеральных Вод занимает особое место среди курортных регионов России благодаря богатству, разнообразию, количеству и ценности минеральных вод, ландшафтных и климатических условий, лечебных грязей. В последнее время возросли темпы освоения ресурсов минеральных вод для целей курорта и промышленного розлива.
Увеличение количества предприятий по розливу минеральной воды и организаций санаторно-курортного типа оказывает существенное влияние на рост объемов водозабора минеральных вод. Нерациональный водоотбор приводит к ухудшению качества подземных вод, изменению их химического состава, температуры. Рост депрессионной воронки может привести к обрушению кровли пласта и исчезновению многих источников. Это касается всех вод, расположенных в регионе Кавказских Минеральных Вод. Существует потенциальная опасность деградации рассматриваемых месторождений минеральных вод.
Важной задачей исследователей является составление прогнозных моделей развития гидролитосфер-ных процессов региона при изменении объемов водозабора на различных участках месторождения. Для этого осуществляется анализ аэрофотоснимков, полученных с помощью беспилотных летательных аппаратов за счет простых линейных алгоритмов.
Для этих целей предлагается использовать технологию Nvidia CUDA с адаптацией под нее математического аппарата анализа фотоснимков. Исходные данные для обработки были получены путем аэрофотосъемки средствами дистанционного зондирования местности беспилотными летательными аппаратами предприятия добычи минеральной воды ОАО «Нарзан» в городе Кисловодске. Программные алгоритмы, приведенные в статье, имеют авторские свидетельства, выданные Федеральным институтом промышленной собственности РФ.
Ключевые слова: функция Грина, тепловое поле, шаг дискретизации, объект управления, анализ, синтез.
Как цитировать эту статью: Первухин Д.А. Параллельный анализ геоданных гидролитосферных пластов минеральной воды Кисловодского месторождения нарзана / Д.А.Первухин, Ю.В.Ильюшин // Записки Горного института. 2016. Т.221. С.706-711. DOI 10.18454/PMI.2016.5.706
Введение. В процессе добычи минеральной воды региона Кавказских Минеральных Вод одной из главнейших задач является нанесение минимального ущерба экологии региона. Это задача, прежде всего, решается с помощью методов очистки и принципиально новых методов добычи минеральной воды. Однако с развитием методов добычи минеральной воды все более актуальным становится вопрос нахождения оптимальных способов поиска месторождений минеральной воды. До недавнего времени в регионе Кавказских Минеральных Вод подобные задачи решались геологическими партиями, проводившими исследования на территории северной части Главного Кавказского хребта. Впоследствии поиск минеральной воды осуществлялся с помощью анализа химического состава рек, берущих свое начало у подножия горы Эльбрус. Этот метод имеет свои недостатки, что связано с тем, что реки образуются за счет огромного количества источников и притоков, впадающих в них (например, река Кич-малка имеет не менее 96 притоков), что меняет химический состав воды. Это существенно затрудняет возможности однозначно идентифицировать класс нарзана (минеральные воды Кисловодского месторождения) и места его залегания.
Для анализа химического состава и температуры воды предложено использовать беспилотные летательные аппараты (БПЛА), которые с помощью аэроспектральной фотографии проводят фотосъемку интересующей местности [1-4]. В последующем данные фотографии проходят обработку. Применение линейного метода обработки фотографий осложняется отсутствием возможности обработки большого объема данных на современных персональных компьютерах. Данную задачу могут решить многопроцессорные и многоядерные системы, однако для обработки данных в многопроцессорных системах персональных компьютеров необходимо использование дополнительных плат расширения, которые существенно увеличивают стоимость системы в целом. С целью снижения стоимости предлагается применение технологии Nvidia CUDA (Compute Unified Device Architecture) для обработки геоданных БПЛА.
Методика исследований. Для обработки данных с БПЛА «Trimble» существует программное обеспечение - фотограмметрический модуль, например Trimble Business Center (TBC) [3]. Для увеличения производительности программного обеспечения необходимо использовать параллельные процессоры.
Вместе с тем эффективное использование современных технических возможностей параллельных процессоров требует принципиально новых подходов к организации вычислений и программи-
ё Д.А.Первухин, Ю.В.Ильюшин
Параллельный анализ геоданных гидролитосферных пластов.
рованию параллельных распределенных алгоритмов. При этом решение задач моделирования, анализа, синтеза и управления распределенными системами можно реализовать на новом, качественно более высоком уровне, применяя методы распараллеливания линейных алгоритмов. Их сущность заключается в разбиении составляющих множителей и слагаемых на компоненты в виде строк двух- и трехмерного массива данных. Далее строки получившегося многомерного массива направляются для решения в арифметико-логическое устройство многопоточного вычислителя. Количество строк многомерного массива зависит от количества арифметико-логических устройств параллельного процессора.
Анализ глубин водных источников по данным аэрофотосъемки. Рассмотрим конкретный пример. Первый блок данных снят с борта БПЛА «ZALA 421-04ф». Данные для исследования были предоставлены ОАО «Газпром космические системы». Блок состоял из 26 маршрутов. Общее число снимков в блоке составило 595. Для получения снимков использовалась предварительно откалибро-ванная цифровая камера «Canon EOS 500D». Высота полета над местностью составила около 500 м, размер пиксела соответствует приблизительно 8 см на местности. На карте местности были измерены и промаркированы 25 опорных точек. Погрешность определения координат опорных точек не превысила 10 см. Общий перепад высот местности протяженностью около 3 км достигал 70 м.
Блок аэрофотосъемки был обработан в автоматическом режиме по упрощенной схеме, без уравнивания и использования опорных точек. Привязка осуществлялась по центрам проекций, трансформирование снимков проводилось в программе PHOTOMOD «GeoMosaic» без учета рельефа. Последующий контроль полученных ортофотопланов по опорным точкам показал расхождения в опорных точках, превышающие 17 м. Такая невысокая точность ортофотоплана обусловлена как большим перепадом высот, так и неточностью измерений центров проекций в полете БПЛА.
Затем блок был подвергнут строгой фотограмметрической обработке. При уравнивании три точки из измеренных опорных считались контрольными. Среднеквадратическая ошибка уравнивания по опорным точкам составила 15; 16; 12 см, по контрольным точкам - соответственно 23; 29 и 57 см. Расхождения на связующих точках составили 8; 14 и 69 см.
В процессе уравнивания было обнаружено, что координаты центров проекций содержат систематическую ошибку, главная из компонент которой составляет 10,5 м по высоте. Среднеквадратиче-ские ошибки на центрах проекций после вычитания систематической ошибки составили 84; 239 и 75 см. Существенно большая ошибка по направлению полета, скорее всего, связана с неточным определением моментов съемки. Большие ошибки на связующих точках, возможно, связаны с неточной калибровкой камеры и накопленной ошибкой при съемке камерой с щелевым затвором. Наибольшие ошибки на связующих точках наблюдаются на краях и в углах снимков. Дальнейшая обработка блока проводилась по стандартной схеме. Был построен рельеф местности в автоматическом режиме и сделано ортотрансформирование с учетом построенного рельефа.
Для анализа результатов создан массив данных, подлежащих обработке. В этом массиве существует огромное количество столбцов, состоящих из разностных характеристик, сумма которых позволяет определить удаленность объекта исследования.
Приведем алгоритм расчета глубины и химического состава водных источников. Данный алгоритм достаточно просто приводится к параллельному виду. Для этого необходимо просуммировать все элементы столбцов - слоев.
Для демонстрации преимущества применения технологии параллельных вычислений предложим две программы, одна из которых будет складывать элементы в простом линейном алгоритме, другая - в параллельном алгоритме применительно к технологии гибридных вычислений Nvidia CUD A [3, 4]. Код программного линейного алгоритма записывается следующим образом
for j := 1 to n do begin for i := 1 to m do begin s[j] := s[j] + a[ij]; sum := sum + a[i,j] end;
Код программного параллельного алгоритма представлен в виде
_global_void add( int *a, int *b, int *c )
int tid = blockldx.x; if (tid < N)
c[tid] = a[tid] + b[tid]; sum := sum + a[tid]
ё Д.А.Первухин, Ю.В.Ильюшин
Параллельный анализ геоданных гидролитосферных пластов.
Из приведенных программных кодов видно, что у линейного алгоритма присутствуют два идентификатора: i - строки и j - столбца. У параллельного алгоритма такой идентификатор задается в настройках массивов, а идентификатор tid задает номер потока данных. Таким образом, первый приведенный программный код создает единый поток данных, производящих суммирующую операцию. А параллельный алгоритм, в свою очередь, создает группу блоков, равных числу ядер процессора, для которых были созданы параллельные массивы данных. Программный код создания такого рода данных имеет следующий вид [1, 2, 5, 6, 13]:
HANDLE_ERROR( cudaMalloc( (void**)&dev_a, sizeof( int ) ) );
HANDLE_ERROR( cudaMalloc( (void**)&dev_b, sizeof( int ) ) );
HANDLE_ERROR( cudaMalloc( (void**)&dev_c, sizeof( int ) ) );
После создания программного комплекса необходимо копировать данные на вычислительный процессор с помощью программного оператора cudaMemcpyDeviceToHost для организации вычислений по векторным типам данных:
cudaMemcpy( dev_a, &z, sizeof(int),cudaMemcpyHostToDevice );
cudaMemcpy( dev_b, &d, sizeof(int),cudaMemcpyHostToDevice) );
cudaMemcpy( dev_c, &k, sizeof( int ), cudaMemcpyHostToDevice );
Вызов модуля расчета производится по схеме с жесткой привязкой к ядру исполнителя: add<<<1,N>>>( dev_t_mas , dev_z, dev_k1, dev_k2, dev_l, dev_a, dev_t_zad, dev_k, dev_d, dev_s);
Это позволило обеспечить функционирование данного модуля исключительно на первом процессоре вычислителя в неограниченном количестве потоков (основное отличие графического процессора от процессора персонального компьютера состоит в наличии у графического процессора сотен арифметико-логических устройств в отличие от одного у процессора персонального компьютера). Изменяя первый параметр функции, можно производить вычисления на других вычислителях, если таковые имеются в системе. Также комбинация данных параметров позволит создавать программные управления для высокопроизводительных систем на базе вычислителей NVidia tesla. Однако в связи с тем, что вычислители tesla обладают четырьмя встроенными процессорами по 2000 ядер каждый, необходимо производить дополнительный расчет потоковых операций. Тем самым создается возможность параллельной обработки группы данных [8, 11] и достигается рост производительности вычислительных сред.
Результаты анализа фотоснимков Кисловодского месторождения, полученных с БПЛА и предоставленных ОАО «Нарзан», с использованием обоих методов идентичны, однако время, затраченное на расчеты, составило: для линейного алгоритма - 15 мин; для параллельного алгоритма - 58 с. (Программные алгоритмы имеют авторские свидетельства о государственной регистрации программы для ЭВМ № 2013661620, 2014613035, 2014613151, выданные Федеральным институтом промышленной собственности РФ.) Здесь стоит отметить, что при решении задачи дополнительные методики расчета, такие как, например, метод сечения, метод разделения переменных и другие, не применялись.
Анализ температурных полей гидролитосферных пластов. При расчете температурных полей используются дифференциальные уравнения в частных производных второго порядка. Для решения настоящей задачи уравнение представлено следующим образом:
дТ(х, y, z, т) — = a
dt
fd 2Т (х, y, z, т) д 2Т (х, y, z, т) 8 2Т (х, y, z, т)
дх1 + dy2 + 8z2 1 + ,
где Т - температура; т - приращение по времени; t - время; а - коэффициент температуропроводности; Q - мощность внутренних источников тепла.
Рассчитав температурное поле на основе алгоритма параллельного синтеза [2, 4], получим дискретный аналог математической модели температурного поля в виде уравнения
( Т (х -1, у, z, t) - 2Т (х, у, z, t) + Т (х +1, у, z, t) ^ Дх2
Т (х, y, z, t) = a At
T(х, y -1, z, t) - 2T(х, y, z, t) + T(х, y +1, z, t)
+ V +
T (х, y, z -1, t) - 2T (х, y, z, t) + T (х +1, y, z + t, t)
+ Az2
+Q.
ё Д.А.Первухин, Ю.В.Ильюшин
Параллельный анализ геоданных гидролитосферных пластов.
Поскольку особенностью Кисловодского месторождения нарзана является наличие на одной географической точке различных видов нарзана (сульфатный, доломитный, общий), необходимо рассчитывать значения температурных полей для каждого слоя нарзана. С этой целью получим дискретные аналоги для трех видов нарзана, расположенных друг над другом и связанных между собой [1, 2, 5-15]: сульфатный
T (x, y, z, t) = a At
f T (x -1, y, z, t) - 2T (x, y, z, t) + T (x +1, y, z, t) ^ Ax2
+ T (x, y -1, z, t) - 2T (x, y, z, t) + T (x, y +1, z, t) +
+
Ay2
T (x, y, z -1, t) - 2T (x, y, z, t) + T (x +1, y, z +1, t)
Az1
+Q;
доломитныи
T (x, y, z, t) = a At
f T (x -1, y, z, t) - 2T (x, y, z, t) + T (x +1, y, z, t) ^ Ax2
+ T (x, y -1, z, t) - 2T (x, y, z, t) + T (x, y +1, z, t) +
+
Ay2
T (x, y, z -1, t) - 2T (x, y, z, t) + T (x +1, y, z +1, t)
Az2
+Q;
общий
T (x, y, z, t) = a At
f T (x -1, y, z, t) - 2T (x, y, z, t) + T (x +1, y, z, t) ^ Ax2
+ T (x, y -1, z, t) - 2T (x, y, z, t) + T (x, y +1, z, t) +
+
Ay2
T (x, y, z -1, t) - 2T (x, y, z, t) + T (x +1, y, z +1, t)
Az2
+Q.
Далее приведем исходный код программного модуля, осуществляющего процесс моделирования. Стоит обратить внимание на то, что процедура начинается с объявления «global». Ключевое слово _global_в объявлении и определении функции обозначает, что функция выполняется на графической карте и что эта функция вызывается из CPU. При вызове этой функции неизменным атрибутом являются скобки <<< >>>. Компилятор nvcc разделяет исходный код на части, которые выполняются отдельно на CPU и GPU. Код, сгенерированный для таких функций, как kernel, выполняется на GPU:
_global_void Cuda_proccedure(double ****dT)
for (r=1; r<=I-1; r++)
for (fi=1; fi<=30; fi++)
for (x=1; x<=1199; x++)
{ if ((r<5)||(r>8)) {a=22.16;}
else a=0.3769;
dT[r] [fi] [x] [timer]=((T[r-1] [fi] [x] [timer] -
2*T[r] [fi] [x] [timer]+T [r+1] [fi] [x] [timer])/r*r+(T[r] [fi-1][x] [timer]-2*T[r] [fi] [x] [timer]+T[r] [fi+1][x] [timer])/fi*fi+(T[r] [fi] [x-1][timer]-2*T[r] [fi] [x] [timer]+T[r] [fi+1][x+1][timer])/x*x); }
Результаты математического моделирования температурных полей на суперкомпьютере и на обычном персональном компьютере одинаковы до пятого знака после запятой, но время, затраченное на расчет с помощью параллельного алгоритма, значительно меньше.
Результаты и их обсуждение. Использование беспилотных летательных аппаратов в качестве аэрофотосъемочной платформы имеет большое будущее. Применение беспилотных летательных аппаратов совместно с программно-аппаратной платформой Nvidia CUDA позволит существенно увеличить скорость обработки данных, получаемых средствами аэрофотосъемки. Однако не следует забывать, что применение технологии Nvidia CUDA дает рост производительности только в случае
ё Д.А.Первухин, Ю.В.Ильюшин
Параллельный анализ геоданных гидролитосферных пластов.
применения параллельных алгоритмов. Когда аэрофотосъемка ведется на большой площади, программно-аппаратные технологии Nvidia CUDA позволяют существенно сократить время расчета как высотных, так и глубинных характеристик за счет одновременного обхода массива данных с нескольких сторон.
Производительность может увеличиться также за счет построения рельефа местности в двунаправленном режиме. Поскольку аппаратная архитектура Nvidia CUDA несет в себе минимум 96 арифметико-логических устройств, данная система может поддерживать минимум 96 направлений режимов обработки растровых изображений. Такое количество арифметико-логических устройств находится в самых дешевых бюджетных версиях графических ускорителей с поддержкой технологии Nvidia CUDA.
При рассмотрении более расширенных версий данной технологии возможно применение графических ускорителей класса «Tesla». Данные ускорители оборудованы четырьмя процессорами по 2000 арифметико-логических устройств в каждом суммарной производительностью до 1 терафлопса (1 триллион операций в секунду). Применение графического ускорителя на локальном персональном компьютере позволяет построить матрицу высот для карт графического объекта площадью 100 га за 6 с. На обычных персональных компьютерах данная задача выполняется в течение суток.
Выводы
1. Анализ аэроспектральной фотографии с использованием параллельных алгоритмов показал существенное увеличение скорости обработки информации, в том числе при расчете глубин и анализе химического состава рек. Это дает возможность более быстрого составления карт пластов местности.
2. При анализе фотографии по температурному полю установлено, что процесс нагрева протекает неоднородно, в зависимости от высоты над уровнем моря. Это, предположительного, связано с наличием термальных источников, а следовательно, наличием горячих доломитных и сульфатных источников нарзана у подножия Главного Кавказского хребта.
3. Моделирование процессов на гибридном суперкомпьютере потребовало 10 мин. Линейный алгоритм на обычном компьютере моделирует данный процесс в течение 6 ч. Высокая скорость моделирования реализуется, главным образом, за счет 96-ядерного процессора, работающего с более высокой суммарной частотой, чем обычный процессор персонального компьютера [1, 5].
4. Проведенный анализ результатов математического моделирования показывает, что процесс параллельного анализа геоданных гидролитосферных пластов минеральной воды можно значительно ускорить за счет применения гибридных суперкомпьютеров. Их использование позволит существенно уменьшить время моделирования переходных характеристик.
ЛИТЕРАТУРА
1. Ильюшин Ю.В. Метод управления температурным полем на основе функции Грина / Ю.В.Ильюшин, И.М.Першин // Записки Горного института. 2015. Т.214. С.57-70.
2. Ильюшин Ю.В. Многопоточный анализ данных сейсморазведки на гибридных суперкомпьютерах / Ю.В.Ильюшин, В.Е.Трушников // Горный информационно-аналитический бюллетень. 2015. № 5. С.230-235.
3. Ильюшин Ю.В. Разработка адаптивной системы управления с распределенным пи-регулятором / Ю.В.Ильюшин, В.Е.Трушников, А.Л.Ляшенко // Горный информационно-аналитический бюллетень. 2014. № 1. С.341-346.
4. Ильюшин Ю.В. Разработка импульсного управления температурным полем буровых шнеков добычи горячей минеральной воды Кисловодского месторождения нарзана / Ю.В.Ильюшин, В.Е.Трушников // Горный информационно-аналитический бюллетень. 2014. № 2. С.172-181.
5. Designing of Distributed Control System with Pulse Control / Yllyushin, D.Pervukhin, O.Afanasieva, A.Klavdiev, S.Kolesnichenko // Middle-East Journal of Scientific Research. 2014. N 21 (3). P.436-439. http://dx.doi.org/10.5829/ idosi.mejsr.2014.21.03.21433
6. Ilyushin Y. Designing of temperature field control system of tunnel kilns of conveyor type // Scientific-technical news of S.Pt.SPI. 2011. N 3 (126). P.67-72.
7. Kolesnikov A. Discharge of a Copper-Magnesium Galvanic Cell in the Presence of a Weak Electromagnetic Field / A.Kolesnikov, Ya.Zarembo, V.Zarembo // Russian Journal of Physical Chemistry. 2007. N 81 (7). P. 1178-1180. http://dx.doi.org/10.1134/s003602440707031x
8. Kolesnikov A. Nonlinear Oscillations Control. Energy Invariants // Journal of Computer and Systems Sciences International. 2009. N 48 (2). P.185-198. http://dx.doi.org/10.1134/S1064230709020038
9. Pleshivtseva Y. The Successive Parameterization Method of Control Actions in Boundary Value Optimal Control Problems For Distributed Parameter Systems / Y.Pleshivtseva, E.Rapoport // Journal of Computer and Systems Sciences International. 2009. N 48 (3). P.351-362. http://dx.doi.org/10.1134/S1064230709030034
10. Rapoport E. Structural Parametric Synthesis of Automatic Control Systems with Distributed Parameters // Journal of Computer and Systems Sciences International. 2006. N 45 (4). P.553-566. http://dx.doi.org/10.1134/S1064230706040071
ё Д.А.Первухин, Ю.В.Ильюшин
Параллельный анализ геоданных гидролитосферных пластов.
11. Solution of Problem of Heating Elements' Location of Distributed Control Objects / Y.Ilyushin, D.Pervukhin, O.Afanasieva, A.Klavdiev, S.Kolesnichenko // Global Journal of Pure and Applied Mathematics. 2016. Vol.12. N 1. P.585-602.
12. Structuring of Inorganic Materials in Weak Rf Electromagnetic Fields / VZarembo, O.Kiseleva, A.Kolesnikov, N.Burnos, K.Suvorov // Inorganic Materials. 2004. N 40 (1). P.86-91. http://dx.doi.org/10.1023/B:INMA.0000012184.66606.59
13. The methods of the synthesis of the nonlinear regulators for the distributed one-dimension control objects / Y.Ilyushin, D.Pervukhin, O.Afanasieva, A.Klavdiev, S.Kolesnichenko // Modern Applied Science. 2015. Vol.9. N 2. P.42-61.
14. Zarembo V. Background Resonant Acoustic Control of Heterophase Processes / V.Zarembo, A.Kolesnikov // Theoretical Foundations of Chemical Engineering. 2006. N 40 (5). P.483-495. http://dx.doi.org/10.1134/s0040579506050058
15. Zinc Electrochemical Reduction on a Steel Cathode in a Weak Electromagnetic Field / A.Kolesnikov, Ya.Zarembo, L.Puchkov, V.Zarembo // Russian Journal of Physical Chemistry. 2007. N 81 (10). P.1715-1717. http://dx.doi.org/10.1134/ s0036024407100330
Авторы: Д.А.Первухин, д-р техн. наук, профессор, [email protected] (Санкт-Петербургский горный университет, Россия), Ю.В.Ильюшин, канд. техн. наук, доцент, [email protected] (Санкт-Петербургский горный университет, Россия).
Статья принята к публикации 02.06.2016.