3

我正在尝试sgemm在 python 中使用 ctypes 来利用 BLAS 中的函数。尝试解决C = A x B以下代码可以正常工作:

no_trans = c_char("n")
m = c_int(number_of_rows_of_A)
n = c_int(number_of_columns_of_B)
k = c_int(number_of_columns_of_A)
one = c_float(1.0)
zero = c_float(0.0)

blaslib.sgemm_(byref(no_trans), byref(no_trans), byref(m), byref(n), byref(k),
               byref(one), A, byref(m), B, byref(k), byref(zero), C, byref(m))

现在我想解这个方程:C = A' x A其中A'是 A 的转置,下面的代码运行没有异常,但返回的结果是错误的:

trans = c_char("t")
no_trans = c_char("n")
m = c_int(number_of_rows_of_A)
n = c_int(number_of_columns_of_A)
one = c_float(1.0)
zero = c_float(0.0)

blaslib.sgemm_(byref(trans), byref(no_trans), byref(n), byref(n), byref(m),
               byref(one), A, byref(m), A, byref(m), byref(zero), C, byref(n))

对于测试,我插入了一个矩阵A = [1 2; 3 4]。正确的结果是C = [10 14; 14 20]sgemm例程会吐出C = [5 11; 11 25]

据我了解,矩阵A不必由我转置,因为算法会处理它。在第二种情况下我的参数传递有什么问题?

任何帮助,链接,文章,建议表示赞赏!

4

2 回答 2

6

Blas 通常使用列主矩阵(如 Fortran),因此A = [1 2; 3 4]意味着

    |1 3|   
A = |   |
    |2 4|

并且结果是正确的(假设您的 Python 库也是如此)。看这个自述文件

于 2011-04-18T16:19:02.190 回答
1

您所指示的结果已sgemm按照您的要求计算了 A*A' 而不是 A'*A 。简单的解决方案是将两个输入切换到函数。

于 2011-04-18T16:10:19.003 回答