LAPACK QR 分解

Question

我正在尝试从 LAPACK zgeqrf 和 zungqr 例程中获取 Q 矩阵。

我有一个 Nw×Nw 复矩阵，其列上有非正交向量。这是我的 C++ 代码。（矩阵命名为vr_tr）

//QR Fact.
complex<double> TAU[Nw*Nw];
zgeqrf_(&Nw, &Nw, vr_tr, &Nw, TAU, &wkopt, &lwork, &info);
lwork = (int)wkopt.real();
work = new complex<double> [lwork];
zgeqrf_(&Nw, &Nw, vr_tr, &Nw, TAU, work, &lwork, &info);
zungqr_(&Nw, &Nw, &Nw, vr_tr, &Nw, TAU, &wkopt, &lwork, &info);
lwork = (int)wkopt.real();
work = new complex<double> [lwork];
zungqr_(&Nw, &Nw, &Nw, vr_tr, &Nw, TAU, work, &lwork, &info);

//Checking if vr_tr * vr_tr' = I
complex<double> one = 1,zer=0,res[Nw*Nw]; char Tchar = 'T', Nchar = 'N';
zgemm_( &Nchar, &Tchar, &Nw, &Nw, &Nw, &one, vr_tr, &Nw, vr_tr, &Nw, &zer, res, &Nw );
for(i=0;i<Nw;i++){
    for(j=0;j<Nw;j++){
        cout<<res[i*Nw+j]<<"\t";
    }
    cout<<"\n\n";
}

我在运行之后得到的不是单位矩阵，因为它应该从 zgeqrf 计算 QR 分解并从 [=12= 获得 Q 矩阵] 和 Q 矩阵必须正交，所以 Q*Q'=I.

这段代码有什么问题？

Answer 1

我用过 zunmqr instead of znugqr，它已经修复了。

虽然老实说我不知道为什么 zungqr 不起作用。

LAPACK QR 分解

LAPACK QR factorization

c++

matrix

lapack