пакетное 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< t> csrlsvlu () для одного b_i, есть ли стандартный способ пакетного выполнения этой операции для нескольких b_i?

Наконец, поскольку у меня нет интуиции для этого, следует ли ожидать ускорения на графическом процессоре для любого из этих параметров, учитывая необходимые накладные расходы? х имеет длину ~ 10000-100000. Благодарю.

6

Решение

Я сам сейчас работаю над чем-то похожим. Я решил в основном обернуть сопряженные градиенты и неполные предварительно подготовленные конъюгаты градиентных решателей уровня 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) Полученное решение ортогонально пустому пространству.

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

Наконец, если ваша матрица имеет Структура циркулянт Вы можете рассмотреть возможность использования решателя на основе быстрого преобразования Фурье (FFT). Числовые решатели на основе БПФ часто дают превосходную производительность в тех случаях, когда они применимы.

1

Другие решения

Если вы не возражаете против использования библиотеки с открытым исходным кодом, вы также можете проверить CUSP:
CUSP Страница быстрого старта

Он имеет довольно приличный набор решателей, включая несколько предобусловленных методов:
Примеры прекондиционера CUSP

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

0