LAPACKE_zheevx() 未能收敛 -- 如何在 C++ 中增加 2*DLAMCH('S') 的 ABSTOL?

LAPACKE_zheevx() failed to converge -- how to increase ABSTOL with 2*DLAMCH('S') in C++?

本文关键字:增加 DLAMCH ABSTOL zheevx LAPACKE C++      更新时间:2023-10-16

这是一个关于在c++中使用LAPACKE_zheevx()函数设置特征值计算收敛的适当公差("abstol")的问题。

当LAPACKE_zheev()在使用默认值"abstol"(即abstol=-1)计算特征值/特征向量时不能收敛时,LAPACK手册说设置abstol=2*DLAMCH('S')。然而,DLAMCH是Fortran函数,我使用的c++不承认它是一个有效的c++函数。谁能帮我如何正确设置"abstol=2*DLAMCH('S')"使用LAPACK与c++(即使用LAPACKE时)?

非常感谢提前!!

背景:LAPACKE是LAPACK(用于数值代数的Fortran库)的c++接口。LAPACKE_zheevx()是LAPACKE的c++接口,用于LAPACK的ZHEEVX()函数。

关键词:LAPACK, LAPACKE, c++, ABSTOL, DLAMCH,收敛性,特征值,特征向量

我对这些库一无所知,但是google一下发现有一个相应的LAPACKE_dlamch()函数:

double LAPACKE_dlamch(char cmach)

因此,您应该能够将LAPACKE_dlamch('S')作为LAPACKE_zheevx()abstol(第12)参数传递:

lapack_int LAPACKE_zheevx (
    int matrix_layout,
    char jobz,
    char range,
    char uplo,
    lapack_int n,
    lapack_complex_double *a,
    lapack_int lda,
    double vl,
    double vu,
    lapack_int il,
    lapack_int iu,
    double abstol,
//  ^^^^^^^^^^^^^
    lapack_int *m,
    double *w,
    lapack_complex_double *z,
    lapack_int ldz,
    lapack_int *ifail
)

除了Jon Purdy的答案(这实际上是正确的):

下面是Fortran中dlamch.f的定义。遍历该函数,并从这里获得Lapack-intrinsic函数,可以看到对于输入's',转换为c++的函数变为:

auto dlamch_s()
{
    auto one = 1.0;
    auto rnd = one;
    auto eps = (one == rnd ? 0.5 : 1.0) * std::numeric_limits<double>::epsilon();
    auto sfmin = std::numeric_limits<double>::min();
    auto small = one / std::numeric_limits<double>::max();
    if(small>=sfmin)
    {
        sfmin = small*(one + eps);
    }
    return sfmin;
}

但是,不要问我为什么需要one == rnd步骤。