Fortran 子例程在 C++ 程序中调用时提供错误结果答案

【问题标题】：Fortran subroutine delivers wrong result when called in C++ programFortran 子例程在 C++ 程序中调用时提供错误结果
【发布时间】：2016-07-22 12:48:59
【问题描述】：

我必须编写一个返回逆矩阵的 Fortran 例程。如果我在 Fortran 程序中运行下面的代码，则逆矩阵是正确的，但是当我从 C++ 代码运行子例程时，我的第一个值是错误的值。这似乎是数据类型或内存的问题。

我做错了什么？

这是子程序：

    subroutine get_inverse_matrix( matrix, rows_matrix, cols_matrix, tmpMatrix, rows_tmpMatrix, cols_tmpMatrix) bind(c)
        use iso_c_binding
        integer :: m, n, lda, lwork, info, size_m
        integer(c_int) :: rows_matrix, cols_matrix, rows_tmpMatrix, cols_tmpMatrix
        real(c_double) :: matrix(rows_matrix, cols_matrix), tmpMatrix(rows_tmpMatrix, cols_tmpMatrix)
        integer, dimension( rows_matrix ) :: ipiv
        real, dimension( rows_matrix )  :: work
        size_m = rows_matrix
        m = size_m
        n = size_m
        lda = size_m
        lwork = size_m
        write(*,*) "Matrix: ", matrix
        !tmpMatrix = matrix
        write(*,*) "Temp matrix: ", tmpMatrix
        ! LU-Faktorisierung (Dreieckszerlegung) der Matrix
        call sgetrf( m, n, tmpMatrix, lda, ipiv, info )
        write(*,*) info
        ! Inverse der LU-faktorisierten Matrix
        call sgetri( n, tmpMatrix, lda, ipiv, work, lwork, info )
        write(*,*) info
        select case(info)
            case(0)
            write(*,*) "SUCCESS"
            case(:-1)
            write(*,*) "ILLEGAL VALUE"
            case(1:)
            write(*,*) "SINGULAR MATRIX"
        end select
    end subroutine get_inverse_matrix

这是 C++ 代码中的声明：

extern "C"
{
void get_inverse_matrix( double *matrix, int *rows_matrix, int *cols_matrix, double *tmpMatrix, int *rows_tmpMatrix, int *cols_tmpMatrix);}

这是来自我的 C++ 程序的调用：

get_inverse_matrix(&lhs[0], &sz, &sz, &res[0], &sz, &sz);

我的程序只使用 3x3 矩阵。如果我通过单位矩阵，结果如下所示：

5.29981e-315 0 0 
0 1 0 
0 0 1

【问题讨论】：

您将声明为c_double 类型的数组传递给期望单精度的lapack 例程，这有时会导致问题。您可以尝试将sgetrf 和sgetri 替换为dgetrf 和dgetri，这有帮助吗？在这种情况下，我也不希望 fortran 能够工作。请注意，您可能还应该在 C++ 代码中显示声明。
我在 C++ 代码中添加了声明。我也按照你的要求替换了，但结果是一个只有零的矩阵，不管我尝试了什么值。
最好使用现代 Fortran90+ 接口模块来 LAPACK，或者至少有一个接口块来检查调用是否有错误。
同样在你的代码中，只要你传入矩阵，但不要将它用于任何事情（除了打印），因为你已经注释掉了 tmpMatrix = matrix 行。这是正确的还是只是一个错字？ matrix 和/或 tmpMatrix（即 lhs 和/或 res）是否已在 C++ 代码中初始化？
感谢您注意到这一点。在尝试修复时，我注释掉了该行。使用线路并且 dgetrf 和 dgetri 调用它有效。非常感谢！

标签： fortran inverse

【解决方案1】：

您将数组声明为类型为real 的类型为c_double，但您使用的是期望单精度输入的lapack 例程（例如c_float）。要解决此问题，您应该将调用 sgetrf 和 sgetri 替换为 dgetrf 和 dgetri。

正如 Vladimir F 在 cmets 中指出的那样，如果您提供接口，这些问题会更容易被发现。

【讨论】：