Пакетное CUDA-решение разреженного Ax=b для различных b
У меня есть разреженная матрица A, и я бы хотел (направить) решение Ax=b. У меня есть около 500 векторов b, поэтому я хотел бы найти соответствующие 500 х. Я новичок в CUDA, поэтому я немного смущен тем, какие варианты у меня есть.
cuSOLVER имеет пакетный прямой решатель cuSolverSP для разреженных A_i x_i = b_i, используя QR здесь. (С LU у меня тоже все будет в порядке, так как A прилично обусловлен.) Однако, насколько я могу судить, я не могу использовать тот факт, что все мои A_i одинаковы.
Будет ли альтернативным вариантом сначала определить разреженную факторизацию LU (QR) на процессоре или графическом процессоре, а затем выполнить параллельную обратную подстановку (соответственно backsub и matrix mult) на графическом процессоре? Если cusolverSp
Наконец, поскольку у меня нет интуиции для этого, следует ли ожидать ускорения на графическом процессоре для любого из этих параметров, учитывая необходимые накладные расходы? х имеет длину ~10000-100000. Благодарю.
2 ответа
Я сам сейчас работаю над чем-то похожим. Я решил в основном обернуть сопряженные градиенты и неполные предварительно подготовленные конъюгаты градиентных решателей уровня 0, которые поставлялись с CUDA SDK, в небольшой класс.
Вы можете найти их в вашем каталоге CUDA_HOME по пути:samples/7_CUDALibraries/conjugateGradient а также /Developer/NVIDIA/CUDA-samples/7_CUDALibraries/conjugateGradientPrecond
По сути, вы должны загрузить матрицу в память устройства один раз (и для ICCG вычислить соответствующий анализ кондиционера / матрицы), а затем вызвать ядро решения с различными векторами b.
Я не знаю, как вы ожидаете, что ваша матричная структура полосы будет выглядеть, но если она симметрична и либо доминирует по диагонали (вне диагональных полос вдоль каждой строки и столбца противоположный знак диагонали, а их сумма меньше, чем у диагонального входа) или положительно определен (нет собственных векторов с собственным значением 0.) тогда CG и ICCG должны быть полезны. С другой стороны, различные многосеточные алгоритмы являются еще одним вариантом, если вы хотите пройти через их кодирование.
Если ваша матрица является только положительно-полуопределенной (например, имеет по крайней мере один собственный вектор с собственным значением, равным нулю), вы все равно можете отказаться от использования CG или ICCG при условии, что: 1) Правая сторона (векторы b) ортогональны нулевому пространству (нулевое пространство означает собственные векторы с собственным значением нуля). 2) Полученное решение ортогонально нулевому пространству.
Интересно отметить, что если у вас есть нетривиальное пустое пространство, то разные числовые решатели могут дать вам разные ответы для одной и той же точной системы. Решения в конечном итоге будут различаться линейной комбинацией нулевого пространства... Эта проблема вызвала у меня много человеко-часов отладки и разочарования, прежде чем я наконец понял, так что хорошо знать об этом.
Наконец, если ваша матрица имеет структуру Circulant Band, вы можете рассмотреть возможность использования решателя на основе быстрого преобразования Фурье (FFT). Числовые решатели на основе БПФ часто дают превосходную производительность в тех случаях, когда они применимы.
есть ли стандартный способ пакетного выполнения этой операции для нескольких b_i?
Один из вариантов - использовать модуль пакетной рефакторизации в CUDA cuSOLVER, но я не уверен, является ли он стандартным.
Модуль пакетной рефакторизации в cuSOLVER обеспечивает эффективный метод решения пакетов линейных систем с фиксированной разреженной матрицей с левой стороны (или матриц с фиксированным шаблоном разреженности, но с переменными коэффициентами) и изменяющимися правыми частями на основе разложения LU. В официальной документации (начиная с CUDA 10.1) можно найти только некоторые частично заполненные фрагменты кода, относящиеся к нему. Полный пример можно найти здесь.
Если вы не возражаете против использования библиотеки с открытым исходным кодом, вы также можете проверить CUSP: CUSP Quick Start Page
Он имеет довольно приличный набор решателей, в том числе несколько предобусловленных методов: Примеры предобусловливателей CUSP
Сглаженное предобусловливание агрегации (вариант алгебраической многосетки), кажется, работает очень хорошо, если у вашего GPU достаточно встроенной памяти для этого.