Armadillo: eigs_gen для наименьшего собственного значения

Я использую eigs_gen броненосца, чтобы найти наименьшее алгебраическое собственное значение разреженной матрицы.

Если я запрашиваю функцию только для наименьшего собственного значения, результат неверен, но если я запрашиваю ее для 2 наименьших собственных значений, результат будет правильным. Код является:

#include <iostream>
#include <armadillo>

using namespace std;
using namespace arma;


int
main(int argc, char** argv)
  {
  cout << "Armadillo version: " << arma_version::as_string() << endl;

  sp_mat A(5,5);

  A(1,2) = -1;
  A(2,1) = -1;
  A(3,4) = -1;
  A(4,3) = -1;

  cx_vec eigval;
  cx_mat eigvec;

  eigs_gen(eigval, eigvec, A, 1, "sr");  // find smallest eigenvalue ---> INCORRECT RESULTS
  eigval.print("Smallest real eigval:");

  eigs_gen(eigval, eigvec, A, 2, "sr");  // find 2 smallest eigenvalues ---> ALMOST CORRECT RESULTS 
  eigval.print("Two smallest real eigvals:");

  return 0;
  }

Моя команда компиляции:

 g++ file.cpp -o file.exe -O2 -I/path-to-armadillo/armadillo-4.600.3/include -DARMA_DONT_USE_WRAPPER -lblas -llapack -larpack

Выход:

Armadillo version: 4.600.3 (Off The Reservation)
Smallest real eigval:
    (+1.000e+00,+0.000e+00)
Two smallest real eigvals:
    (-1.000e+00,+0.000e+00)
    (-1.164e-17,+0.000e+00)

Любая идея о том, почему это происходит и как это преодолеть, приветствуется.

Примечание: второй результат является почти почти правильным, потому что мы ожидаем -1, -1 в качестве двух самых низких собственных значений, но, возможно, повторные собственные значения игнорируются.

Обновление: включая построение тестовой матрицы, которая после изменений Райана, чтобы включить опцию "sa" в библиотеку, кажется, не сходится:

    #define ARMA_64BIT_WORD
    #include <armadillo>
    #include <iostream>
    #include <vector>
    #include <stdio.h>

    using namespace arma;
    using namespace std;

    int main(){

    size_t l(3), ls(l*l*l);
    sp_mat A = sprandn<sp_mat>(ls, ls, 0.01);
    sp_mat B = A.t()*A;

    vec eigval;
    mat eigvec;
    eigs_sym(eigval, eigvec, B, 1, "sa");

    return 0;

    }

Интересующие нас размеры матрицы намного больше, например ls = 8000 - 27000и не совсем матрица, построенная здесь, но я предполагаю, что проблема должна быть такой же.

Источник

2 ответа

Решение

Я считаю, что проблема в том, что вы работаете eigs_gen() (который вызывает DNAUPD) на симметричной матрице. ARPACK отмечает, что DNAUPD не предназначен для симметричных матриц, но не определяет, что произойдет, если вы все равно будете использовать симметричные матрицы:

ПРИМЕЧАНИЕ. Если линейный оператор "OP" является действительным и симметричным относительно действительной положительной полуопределенной симметричной матрицы B, т. Е. B*OP = (OP')*B, то вместо этого следует использовать подпрограмму ssaupd.

(с http://www.mathkeisan.com/usersguide/man/dnaupd.html)

Я изменил внутренний код Armadillo, чтобы передать "sa" (наименьший алгебраический) вызовам ARPACK в eigs_sym() (sp_auxlib_meat.hpp), и я смог получить правильные собственные значения. Я отправил патч для апстрима, чтобы сделать поддержку "sa" и "la" доступной для eigs_sym(), которая, я думаю, должна решить вашу проблему, как только выйдет новая версия (или в какой-то момент в будущем).

Источник

Проблема в повторяющихся собственных значениях; если я изменю первые два элемента матрицы на

A(1,2) = -1.00000001;
A(2,1) = -1.00000001;

ожидаемые результаты получены.

Источник
Другие вопросы по тегам