使用GSL进行矩阵乘法

Question

我有两个GSL矩阵，AT和A:

gsl_matrix *A;    /* coefficient matrix A     */
gsl_matrix *AT;   /* coefficient matrix A'    */

AT   = gsl_matrix_alloc(nc, nr); /* Data matrix */
A    = gsl_matrix_alloc(nr, nc); /* Data matrix */

/* Initialize A */
for(i = 0; i < nr; i++){
  gsl_matrix_set(A, i, 0, 1.0); 
}

for(i = 0; i < nr; i++){
  for(j = 1; j < nc; j++){
    gsl_matrix_set(A, i, j, 1.0 / (double)(i + j + 1)); 
  }
}

gsl_matrix_transpose_memcpy(AT, A); 

我希望将这些矩阵相乘并将结果存储在矩阵ATA中，但我无法理解GSL BLAS文档。

我已经初始化了ATA矩阵：

gsl_matrix *ATA;  /* coefficient matrix A'A   */
ATA  = gsl_matrix_alloc(nc, nc); /* Data matrix */  

我发现我可以使用来乘以复杂矩阵，但这些矩阵并不是复杂的。那么我该如何处理呢？

更新：

我已经尝试过：

gsl_blas_dgemm(CblasNoTrans, CblasNoTrans, 1.0, A, AT, 0.0, ATA); 

它导致了以下错误：

gsl: blas.c:1354: 错误：无效长度 默认的GSL错误处理程序被调用。 异常终止（核心已转储）

- Sara Fuerst

Answer 1

我把矩阵乘法反过来了。正确的做法是：

gsl_blas_dgemm(CblasNoTrans, CblasNoTrans, 1.0, AT, A, 0.0, ATA); 

- Sara Fuerst