LAPACKE C ++实矩阵求逆

Question

我正在尝试将C ++ LAPACKE中的真实矩阵求逆。对于复杂矩阵，我具有相同的功能，并且可以正常工作。但是实际情况给出了错误的答案。这是我的功能：

void inv(std::vector<std::vector<double>> &ans, std::vector<std::vector<double>> MAT){

    int N = MAT.size();

    int *IPIV = new int[N];

    double * arr = new double[N*N];
    for (int i = 0; i<N; i++){
        for (int j = 0; j<N; j++){
            int idx = i*N + j;
            arr[idx] = MAT[i][j];
        }
    }

    LAPACKE_dgetrf(LAPACK_ROW_MAJOR, N, N, arr, N, IPIV);
    LAPACKE_dgetri(LAPACK_ROW_MAJOR, N, arr, N, IPIV);

     for (int i = 0; i<N; i++){
        for (int j = 0; j<N; j++){
            int idx = i*N + j;
            ans[i][j] = arr[idx];
        }
    }

    delete[] IPIV;
    delete[] arr;
}

我尝试将24 x 24的双精度矩阵求反。尽管程序似乎快要存在了，但是逆函数还不完全存在，它与python linalg逆函数给我的有很大不同（python就在这里，因为我将矩阵与逆相乘，结果非常接近于恒等式）。在LAPACKE输出中，我将矩阵乘以其逆，得到对角线为1，但非对角线的值最高为0.17，与0相比，它是巨大的。有没有办法使LAPACKE程序提供更好的结果？谢谢！

Answer 1

对于行列式较大的矩阵，您可以重新缩放输入，计算逆，然后重新缩放输出。这是一个非常简单的Python示例，您的比例因子应为1/25左右，以便获得（1/25）²⁴ = 2.8e-34的总乘数，从而使输入矩阵的行列式约为1000。]

import numpy as np

scale = 0.5

i = np.array([[1,2],[3,4]])
print(i)
print(np.linalg.det(i))
print("-----------------------------------")
x = np.multiply(i, [scale]) # rescale matrix
print(x)
print(np.linalg.det(x))     # determinant should be less
print("-----------------------------------")
y = np.linalg.inv(x)
print(y)
print(np.linalg.det(y))
print("-----------------------------------")
o = np.multiply(y, [scale]) # rescale matrix
print(o)
print(np.linalg.det(o))
print(np.dot(i, o))