使用 gsl、Blas 和 Lapack 计算 (Aᵀ×A)*(Bᵀ×B) 矩阵

Question

让我有 2 个对称矩阵：

A = {{1,2}, {2,3}}
B = {{2,3},{3,4}}

我可以使用 gsl、Blas 和 Lapack 计算矩阵 (A^T×A)*(B^T×B) 吗？

我正在使用

gsl_blas_dsyrk(CblasUpper, CblasTrans, 1.0, A, 0.0, ATA);
gsl_blas_dsyrk(CblasUpper, CblasTrans, 1.0, B, 0.0, BTB);
gsl_blas_dsymm(CblasLeft, CblasUpper, 1.0, ATA, BTB, 0.0, ATABTB); // It doesn't work

它returns：

(Aᵀ·A) = ATA = {{5, 8}, {0, 13}} -- ok, gsl_blas_dsyrk returns symmetric matrix as upper triangular matrix.
(Bᵀ·B) = BTB = {{13, 8}, {0, 25}} -- ok.
(Aᵀ·A)·(Bᵀ·B) = ATABTB = {{65, 290}, {104, 469}} -- it's wrong.

Answer 1

对称BTB问题就解决了

如您所见，对称矩阵的上三角部分由 dsyrk() 计算。然后应用 dsymm()。根据dsymm(), the following operation is performed since the flag CblasLeft的定义使用：

C := alpha*A*B + beta*C

where alpha and beta are scalars, A is a symmetric matrix and B and C are m by n matrices.

的确，B矩阵是一般矩阵，不一定是对称矩阵。结果，ATA 乘以 BTB 的上三角部分，因为未计算 BTB 的下三角部分。

将 BTB 对称化，问题将得到解决。为此，for 循环是一个直接的解决方案，请参阅

使用 gsl、Blas 和 Lapack 计算 (Aᵀ×A)*(Bᵀ×B) 矩阵

Compute (Aᵀ×A)*(Bᵀ×B) matrix using gsl, Blas and Lapack

matrix

blas

lapack

matrix-multiplication

gsl