Инструмент для диагонализации больших матриц

Question

Я хочу вычислить диффузионное ядро, которое включает в себя выражение exp (b * A), где A - большая матрица. Чтобы поиграть со значениями b, я бы хотел диагонализовать A (чтобы exp (A) выполнялась быстро).

Моя матрица размером около 25k x 25k, но она очень разреженная - только около 60k значений не равны нулю. Функция Matlab "eigs" запускает нехватку памяти, как и октава "eig" и R "eigen". Есть ли инструмент для поиска разложения больших разреженных матриц?

Не знаю, если это уместно, но A - это матрица смежности, поэтому она симметрична и имеет полный ранг.

Answer 1

Вы пробовали SVD, svds для разреженной матрицы в Matlab.

РЕДАКТИРОВАТЬ: еще одна вещь, не делайте SVD полного ранга, так как размерность велика, используйте маленький ранг, скажем 500, чтобы ваше решение поместилось в памяти. Это отсекает маленькие собственные значения и их векторы. Таким образом, это не повредит вашей точности.

Answer 2

Если у вас есть доступ к 64-битной машине и октаве, скомпилированной с поддержкой 64-битной, вы можете обойти эту проблему.

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

Например, вы можете запустить

ulimit -u unlimited

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

ulimit -s [stacksize]

, чтобы увеличить ограничение размера стека.

Answer 3

Рассматривали ли вы следующее свойство : exp (A * t) = L ^ (- 1) {(sI-A) ^ (- 1)} где L ^ (- 1) обратное преобразование Лапласа? - при условии, что вы можете инвертировать (sI-A)

Answer 4

В R вы можете проверить пакет igraph и функцию arpack, которая является интерфейсом к библиотеке ARPACK .

Answer 5

В Octave есть splu, которая выполняет разложение lu для разреженных матриц. Я не уверен, может ли он справиться с 25k x 25k, но стоит попробовать.

В качестве альтернативы, если ваша матрица структурирована следующим образом: A = [B нули; нули C], вы можете по диагонали B и C по отдельности и объединить их в одну матрицу. Я предполагаю, что вы можете сделать что-то подобное для EIG.