Гужва А.Г. Использование нейронных сетей для решения обратной задачи Магнитотеллурического Зондирования

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

Задача Магнитотеллурического зондирования (МТЗ) Обратная задача МТЗ ( -> ): A k =, Г к R Nk, R Mk A k – заданный дискретный оператор прямой задачи = ( 1 … Mk ) – вектор характеристик МТ поля, измеренных на поверхности Земли = ( 1 … Nk ) – вектор макропараметров среды Г k – область допустимых значений k – класс разрезов

Задача Магнитотеллурического зондирования (МТЗ) Обратная задача МТЗ ( -> ): A k =, Г к R Nk, R Mk Система нелинейных уравнений относительно Неустойчива и некорректна – Вид уравнений – Размерность данных Чудовищная размерность данных

Задача Магнитотеллурического зондирования (МТЗ) Прямая задача МТЗ ( -> ): A k =, Г к R Nk, R Mk Устойчива и корректна Считается на суперкомпьютерах Эталонная база решений прямой задачи: – Задаем 1, считаем 1 =A k 1, знаем ( 1, 1 ) – … – Задаем N, считаем N =A k N, знаем ( N, N )

Задача Магнитотеллурического зондирования (МТЗ) Приближенное решение обратной задачи – * appr ( 1 … N, a 1,…,a J ) – a i – неопределенные коэффициенты – appr – заданная функция некоторого вида – Решение задачи нелинейной регрессии на эталонной базе решений прямой задачи Полученное решение есть приближенный обратный оператор задачи

Задача Магнитотеллурического зондирования (МТЗ), 2D 336 блоков ( ) Электропроводность <

Задача МТЗ. Пример распределения

Задача нелинейной регрессии Дано: Задан конечный обучающий набор векторов X = {x 1 …x n }, x i R M, i = 1..n Имеется зависимость y = y*(x), y*: R M R P, для которой известны y i = y*(x i ) для x i X Задано разбиение X на тренировочный набор данных X Trn и тестовый набор данных X Tst Задан вид семейства нелинейных параметрических функций M переменных y = f(x, )

Задача нелинейной регрессии Требуется: Аппроксимация y* с помощью f(x, ) Найти решение системы уравнений по y*(x 1 ) = f(x 1, ) … y*(x m ) = f(x m, ) при котором достигался бы минимум функционала ошибки: E(, X Tst ) = x X Tst (y*(x) - f(x, )) 2, x i X Trn

Нейронные сети. Многослойные персептроны МСП - семейство параметрических вектор-функций МСП с 1 скрытым слоем (f: R M R P ) – ( ) С 1 – ограниченная монотонно возрастающая действительная функция – и u – матрицы весов (N x (M + 1) и P x (N + 1)) – k = 1..P – N – характеризует сложность конструкции

Традиционный способ решения задачи нелинейной регрессии для многослойных персептронов (обучение сети) Оптимизационная задача по подбору элементов матриц весов Решение методом наименьших квадратов путем минимизации функционала ошибки: E(, X Trn ) = x X Trn (y*(x) - f(x, )) 2 min E(, X Tst ) используется для контроля останова процесса обучения нейронной сети Например, используется градиентный спуск

Задача МТЗ. 1 большая обратная задача МТЗ – «В лоб» не решается – Размерность входного вектора 6552 – Размерность выходного вектора небольших задач нелинейной регрессии – Размерность входного вектора 1648 – Размерность выходного вектора 1

Задача МТЗ. Формально. Решение ~7000 задач нелинейной регрессии – Вид аппроксимирующей функции y(x, u, ): где u, – матрицы неопределенных коэффициентов Для всех задач нелинейной регрессии входной набор данных одинаков – Размерность вектора входных данных x = 1648 Выходные наборы данных уникальны – Размерность вектора выходных данных y = 1

Задача МТЗ. Формально примеров эталонной базы решений Решение задачи нелинейной регрессии (подбор u, ) – Минимизация функционала ошибки вида где W = {u, } – условное обозначение матриц весов многослойного персептрона y(x, W) Например, использование метода градиентного спуска

CPU. Градиентный спуск Минимизация функционала ошибки вида Размерность W ~ 10 4 Градиентный спуск: 1.W=W 0 2.W=W - *grad E(W), где ~10 -1 – goto 2 Распараллеливается на MPI

CPU 6780 небольших задач нелинейной регрессии – 5 AMD Athlon 64 x2 Dual GHz – 1 ядро того же уровня – Итого 11 ядер Вычисления во float-ах Своя спец. библиотека Итого - 2 месяца вычислений – С перерывами на выходные

CUDA. Распараллеливание Нет – Убыстрение 1 вычисления y(x, u, ), y/ u, y/ Да! – Параллельное вычисление y(x 1, u, ), …, y(x N, u, ) – Параллельное решение нескольких задач нелинейной регрессии y(x, u 1, 1 ), …, y(x, u M, M ) y(x 1, u 1, 1 ), …, y(x N, u M, M )

CUDA. Реализация Функционал ошибки – E(W) = i=1,N (y(x i,W) - y i ) 2 grad E(W) = 2 * i=1,N (y(x i,W) - y i ) * y/ W(x i ) Одновременное решение нескольких задач нелинейной регрессии – grad E 1 (W 1 ) = 2 * i=1,N (y(x i,W 1 ) - y i ) * y/ W |1 (x i ) ……… ……… ………. – grad E k (W k ) = 2 * i=1,N (y(x i,W k ) - y i ) * y/ W |k (x i )

CUDA. Реализация Использование CUBLAS для вычисления значений функции (cublasSgemm) ………… ……… Аналогично с производными

CUDA. Реализация Эффект за счет резкого уменьшения числа операций работы с памятью – Необходимость переноса больших объемов данных между узлами – CUDA: наличие shared memory Основное время при проектировании – Расположение данных в памяти Код спроектирован с нуля Вычисления во float-ах Вид аппроксимационной функции

Результаты CPU: AMD Athlon 64 x2 Dual GHz В GPU-реализации не используется Zero copy memory

GPU 6780 небольших задач нелинейной регрессии – 1 Tesla C2050 – 256 параллельно решаемых задач нелинейной регрессии – Итого: ~13 часов – На разогнанном GTX 580 ~ 6 часов (ориентировочно)

«Грабли». Учитывать при проектировании Доступные ресурсы для ядра – Регистры, shared memory, потоки RTFM: Zero copy, pinned memory, … – А также прочие полезные «фичи» Bottleneck для этой задачи – Пропускная способность внутренней шины данных (~100 GB/s)

«Грабли». Раскрытие циклов //// Ядро с циклом __global__ void kernel(int N) {... for (int i = 0; i < N; i++)//// Цикл...//// Тело цикла... } //// Этот код быстрее! template __global__ void kernel() {... for (int i = 0; i < N; i++) //// Цикл...//// Тело цикла... } N – заранее известно N может принимать огр. число значений

«Грабли». Использование CUDA Primitives CUDA parallel primitives – Reduction – Scan – И т.д. В виде отдельных процедур – CUDPP Внедрить в собственное ядро!

«Грабли». Если глючит, то… nvidia-smi cuda-memcheck Проверять cudaGetLastError(); В код повставлять cudaThreadSynchronize(); Выделение shared памяти __global__ void kernel() { __shared__ float klmn[100]; //// в klmn – не нули!!!!!!... }

«Грабли». Если глючит, то… Воспроизводимость глюка – Запуск процесса в цикле – Каждый раз сравнение результатов, в т.ч. промежуточных массивов данных Компилятор nvcc RTFM, форумы, коллеги

«Грабли». Если глючит, то… Если ничего не помогло, то…

Распараллеливание Вычисление функции y(x, a) Варианты действий: – Распараллеливание y(x, a) – Одновременно y(x 1, a), …, y(x N, a) – Одновременно y(x, a 1 ), …, y(x, a M ) – Одновременно y(x 1, a 1 ), …, y(x N, a M )

Спасибо за внимание!