1

我有代码进行数值优化(大型程序)。在 gdb 中,它以

Program received signal EXC_BAD_ACCESS, Could not access memory.
Reason: 13 at address: 0x0000000000000000
[Switching to process 8673 thread 0x1203]
0x00007fff95d2a9fb in small_malloc_from_free_list ()

在回溯期间,我被告知以下行包含错误,但无法找出问题所在。

177     U  = calloc(cols*cols,sizeof(double));

包含该行的代码是

int fcn_invXX(double *X, int rows, int cols, double *invXX)
{
    //---------------------------------------------------------------
    // Declarations
    //---------------------------------------------------------------

    // Counters and info
    int i, j, info;

    // Upper triangular placeholder, U, and X'X
    double *U, *XX;

    // Indicator for upper triangular
    char s;

    // Alpha and beta for matrix multiplication - see dgemm documentation
    double alpha, beta;


    //---------------------------------------------------------------
    // Memory Allocation
    //---------------------------------------------------------------
    U  = calloc(cols*cols,sizeof(double));
    XX = calloc(cols*cols,sizeof(double));


    //---------------------------------------------------------------
    // Compute inverse of X'X
    //---------------------------------------------------------------
    alpha = 1.0;
    beta  = 0.0;

    // Get X'X
    cblas_dgemm(CblasColMajor, CblasTrans, CblasNoTrans, cols, cols, rows, alpha, X, rows, X, rows, beta, XX, cols);
    // CblasColMajor means that data are ordered in columns first, so if X is 4x4, X[0] and X[1] are first column
    // CblasTrans means that the matrix should be transposed before multiplication. Here first matrix is transposed
    // CblasNoTrans means that matrix should not be transposed. Here second matrix is not transposed
    // X is the first matrix, and here, also the second
    // XX stores the product

    for (i = 0; i < cols*cols; i++)
    {
        U[i] = XX[i];
    }

    // Upper Cholesky factorization
    s = 'U'; // upper
    dpotrf_(&s, &cols, U, &cols, &info);
    // info = LAPACKE_dpotrf(LAPACK_COL_MAJOR, s, cols, U, cols);

    for (i = 0; i < cols*cols; i++)
    {
        invXX[i] = U[i];
    }

     // Get inverse
    dpotri_(&s, &cols, invXX, &cols, &info);
    // info = LAPACKE_dpotri(LAPACK_COL_MAJOR, s, cols, invXX, cols);

    for (i = 0; i < cols; i++)
    {
        for (j = 0; j < cols; j++)
        {
            if (i > j)
            {
                invXX[i + j * cols] = invXX[j + i * cols];
            }
        }
    }


    //---------------------------------------------------------------
    // Free Memory
    //---------------------------------------------------------------
    free(U);
    free(XX);

    return info;

}

希望你能帮忙。

4

0 回答 0