对角化一个2×2的厄米矩阵很简单,它可以用解析的方法来完成。然而,当涉及到计算超过&>106次的特征值和特征向量时,重要的是尽可能高效地进行。特别是如果非对角线元素可以消失,就不可能对特征向量使用一个公式:一个if-statement是必要的,这当然会减慢代码的速度。因此,我认为使用特征值仍然是一个不错的选择,特征值说明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);
}
直接利用厄米矩阵的特征值和特征向量的公式以及一个if语句来检验off对角线是否为零,写出不带本征的循环,比写出循环快5倍。是我没有正确地使用本征,还是这样的开销是正常的?是否有其他的Lib.s是针对小的自伴随矩阵优化的?
默认情况下,使用迭代方法。要对2x2和3x3使用分析版本,必须调用computeDirect函数:
ces.computeDirect(X);
但它不太可能比解析公式的实现更快。