Диагонализация эрмитовой матрицы 2x2 проста, это можно сделать аналитически. Однако, когда дело доходит до вычисления собственных значений и собственных векторов более чем в 10 ^ 6 раз, важно делать это как можно более эффективно. Особенно, если недиагональные элементы могут исчезнуть, невозможно использовать одну формулу для собственных векторов: необходим оператор if, который, конечно, замедляет код. Таким образом, я подумал, что использование Eigen, где сказано, что диагонализация матриц 2x2 и 3x3 оптимизирована, по-прежнему будет хорошим выбором:
с использованием
const std::complex<double> I ( 0.,1. );
inline double block_distr ( double W )
{
return (-W/2. + rand() * W/RAND_MAX);
}
тестовый цикл будет
...
SelfAdjointEigenSolver<Matrix<complex< double >, 2, 2> > ces;
Matrix<complex< double >, 2, 2> X;
for (int i = 0 ; i <iter_MAX; ++i) {
a00=block_distr(100.);
a11=block_distr(100.);
re_a01=block_distr(100.);
im_a01=block_distr(100.);
X(0,0)=a00;
X(1,0)=re_a01-I*im_a01;
//only the lower triangular part is referenced! X(0,1)=0.; <--- not necessary
X(1,1)=a11;
ces.compute(X,ComputeEigenvectors);
}
Написание цикла без Eigen с использованием непосредственно формул для собственных значений и собственных векторов эрмитовой матрицы и оператора if для проверки того, равна ли недиагональ нулю, происходит в 5 раз быстрее. Я неправильно использую Eigen или такие накладные расходы - это нормально? Существуют ли другие библиотеки lib.s, оптимизированные для небольших самосопряженных матриц?
