MKL能求解非对称矩阵部分特征值的问题吗？

木木林

    自己写了个有限元程序，最后需要求解一个大型（大概1000*1000）非对称稀疏矩阵的部分特征值及其特征向量，如最大的10个特征值与对应的特征向量。最近研究了下  Intel MKL  Fortran 的参考手册，发现里面有两个地方提到了特征值的求解：


1.LAPCK---->LAPACK Least Squares and Eigenvalue Problems Routines  

   这里面能够求解四种特征值问题：对称EVPs，对称正定的广义EVPs，非对称EVPs和广义非对称EVPs。还提供了计算Computationl Routines 和Driver Routines 两种版本。应该用的是QR迭代法。

2.Extened Eigensolver Routines
  这里面提到了一种Feast 算法，支持CSR和BSR等压缩矩阵格式，可以求解最大，最小特征值，和在计算给定范围内的特征值和特征向量，但是似乎只支持Symmetric和Hermitian两种性质的矩阵。

那么现在的问题是：
     LAPACK的特征值求解器不支持CSR和BSR等压缩存储格式，也不支持求解部分特征对；而Extended Eigensolver 又不支持非对称矩阵。那么
     MKL能够求解  非对称  稀疏矩阵的   最大的若干个特征值及其特征向量吗？

Jackdaw

?ggev
Computes the generalized eigenvalues, and the left
and/or right generalized eigenvectors for a pair of
nonsymmetric matrices.


可以求所有的特征值和特征向量，1000阶很快

木木林

Jackdaw 发表于 2019-6-5 13:05
?ggev
Computes the generalized eigenvalues, and the left
and/or right generalized eigenvectors for a ...

感谢回复。
我是用1D FEM求一个旋转梁的流固耦合问题，实际中我一般只需要前10-20个特征值和特征向量分别表示固有频率和模态形状够了。

ggev是LAPACK里面的driver routine。这个LAPACK里面的eigensolver 用的是QR迭代法，需要用稠密阵存储矩阵，不像Extended eigensolver 一样支持稀疏矩阵存储格式，如CSR。而且后期有限元网格会加密，可能带来效率问题和特征向量的存储开销问题。如10000*10000，光所有特征向量就需要存储一亿个浮点数，太夸张了。

目前比较清楚的情况是：

   目前的MKL里的Extended eigensolver 用的是FEAST 2.0, 所以不支持非对称和非Hermitian矩阵。
   2015年发布的FEAST 3.0是支持非对称矩阵的，网上看到Intel的人2016年就说注意到了，但直到现在也没加到MKL。

我尝试下看能不能让MKL调FEAST 3.0。

pasuka

木木林 发表于 2019-6-5 23:52
感谢回复。
我是用1D FEM求一个旋转梁的流固耦合问题，实际中我一般只需要前10-20个特征值和特征向量分 ...

直接用ARPACK不是挺好的嘛！
传送门
https://github.com/opencollab/arpack-ng

木木林

pasuka 发表于 2019-6-7 06:58
直接用ARPACK不是挺好的嘛！
传送门
https://github.com/opencollab/arpack-ng

   感谢回复。
   之前对特征值问题完全不了解，又因为用的是windows，只是对MKL熟悉点，所以一直想在MKL的框架里解决问题，同时也是想利用MKL中高度硬件优化的BLAS、SPBLAS和LAPCK库。  

   最近折腾了一段时间后，也了解到了ARPACK，对于非对称、非Hermitian特征值问题总结如下，以供后人学习参考：

    1. 中小型特征值问题
        这类问题可以采用QR迭代法直接求解所有特征值和特征向量，只支持稠密矩阵，主要是LAPACK中的相应subroutine用的就是这类方法。

    2. 大型特征值问题的部分特征值问题。
        2.1 ARPACK
             ARPACK采用Krylov 子空间方法和implicit restart 技术，大概在1996年前后开发完成。
        2.2 Feast
            Feast借鉴量子力学中的围道积分法（contour integration），对由围道（contour）在复平面中所围的区域内进行特征值搜索，可以同时设置多个围道，并行效果较好。Feast 2.0只支持对称和Hermitian问题，Feast 3.0 支持所有矩阵问题（对称，非对称，Hermitian，非Hermitian  , dense  , banded, sparse）。Feas 3.0的开发人员中有intel的人，所以很多BLAS和LAPACK库默认是从MKL中调用的，稀疏特征值问题还用到了MKL中的PARDISO模块，因此编译的时候最好用Intel Fortran 编译器，大概语法为：
    ifort  /Qmkl   /fpp  /O2  /c  /I.  *.f90

具体理论可参考Feast 2.0【1】和3.0【2】的原始论文：
  【1】 Eric Polizzi，2009，A Density Matrix-based Algorithm for Solving Eigenvalue Problems
  【2】JAMES KESTYNy, ERIC POLIZZIz, AND PING TAK PETER TANG，2016，FEAST EIGENSOLVER FOR NON-HERMITIAN PROBLEMS
  

   

  

Jackdaw

木木林 发表于 2019-6-10 20:07
感谢回复。
   之前对特征值问题完全不了解，又因为用的是windows，只是对MKL熟悉点，所以一直想在MKL ...

厉害，我涉及特征值的问题只是为了求谱条件数，用了MKL里GGSV，发现慢，改用了幂法反幂法

pasuka

木木林 发表于 2019-6-10 20:07
感谢回复。
   之前对特征值问题完全不了解，又因为用的是windows，只是对MKL熟悉点，所以一直想在MKL ...

lz的口气好大哦！
1、大型特征值，多大才算大呢？若是大型稠密矩阵呢？
2、Lanczos法、子空间迭代法、Jacob-Davidson法总不能无视吧？
3、FEAST方法的点子很新颖，将子空间迭代法扩展到复数域。可是判断区间内特征根数量有啥好法子嘛？
4、PARDISO是业界标杆，可是SuperLU、UMFPACK、SPQR也蛮好玩的
5、学习有限元开发，为啥不用C++的Deal II

chdlkl

您好，请问压缩存储稀疏对称矩阵后，求矩阵特征值的问题解决了没

asri623

试了一下在用FEAST进行区间特征值求解，发现效率比Lanczos法高，但对在有些问题中经常找不到特征值，搞不懂原因？请指教！