c - 在C中将矩阵相乘

标签 c r matrix

我正在尝试使用C中的循环将几个矩阵相乘。我在R中获得了预期的答案,但在C中无法获得期望的答案。我怀疑问题与+=函数有关,该函数似乎在循环的第一次迭代后将乘积的值加倍。

我对C不太熟悉,也无法用会返回预期答案的+=函数代替。

感谢您的任何建议。

首先,这是返回所需答案的R代码:

B0 = -0.40
B1 =  0.20

mycov1 = exp(B0 + -2 * B1) / (1 + exp(B0 + -2 * B1))
mycov2 = exp(B0 + -1 * B1) / (1 + exp(B0 + -1 * B1))
mycov3 = exp(B0 +  0 * B1) / (1 + exp(B0 +  0 * B1))
mycov4 = exp(B0 +  1 * B1) / (1 + exp(B0 +  1 * B1))

trans1 = matrix(c(1 - 0.25 - mycov1,     mycov1,   0.25 * 0.80,              0,
                                  0,   1 - 0.50,             0,    0.50 * 0.75,
                                  0,          0,             1,              0,
                                  0,          0,             0,              1), 
               nrow=4, ncol=4, byrow=TRUE)

trans2 = matrix(c(1 - 0.25 - mycov2,     mycov2,   0.25 * 0.80,              0,
                                  0,   1 - 0.50,             0,    0.50 * 0.75,
                                  0,          0,             1,              0,
                                  0,          0,             0,              1), 
               nrow=4, ncol=4, byrow=TRUE)

trans3 = matrix(c(1 - 0.25 - mycov3,     mycov3,   0.25 * 0.80,              0,
                                  0,   1 - 0.50,             0,    0.50 * 0.75,
                                  0,          0,             1,              0,
                                  0,          0,             0,              1), 
               nrow=4, ncol=4, byrow=TRUE)

trans4 = matrix(c(1 - 0.25 - mycov4,     mycov4,   0.25 * 0.80,              0,
                                  0,   1 - 0.50,             0,    0.50 * 0.75,
                                  0,          0,             1,              0,
                                  0,          0,             0,              1), 
               nrow=4, ncol=4, byrow=TRUE)

trans2b <- trans1  %*% trans2
trans3b <- trans2b %*% trans3
trans4b <- trans3b %*% trans4
trans4b

#
# This is the expected answer
#
#            [,1]      [,2]      [,3]      [,4]
# [1,] 0.01819965 0.1399834 0.3349504 0.3173467
# [2,] 0.00000000 0.0625000 0.0000000 0.7031250
# [3,] 0.00000000 0.0000000 1.0000000 0.0000000
# [4,] 0.00000000 0.0000000 0.0000000 1.0000000
#


这是我的C代码。 C代码相当长,因为我不太了解C来提高效率:

#include <string.h>
#include <stdio.h>
#include <stdlib.h>
#include <math.h>

char quit;

int main(){

int i, j, k, ii, jj, kk ;

double B0, B1, mycov ;

double trans[4][4]     = {0} ;

double prevtrans[4][4] = {{1,0,0,0},
                          {0,1,0,0},
                          {0,0,1,0},
                          {0,0,0,1}};

B0 = -0.40 ;
B1 =  0.20 ;

for (i=1; i <= 4; i++) {

          mycov = exp(B0 + B1 * (-2+i-1)) / (1 + exp(B0 + B1 * (-2+i-1))) ;

          trans[0][0] =     1 - 0.25 - mycov ;
          trans[0][1] =                mycov ;
          trans[0][2] =          0.25 * 0.80 ;
          trans[0][3] =                    0 ;

          trans[1][0] =                    0 ;
          trans[1][1] =             1 - 0.50 ;
          trans[1][2] =                    0 ;
          trans[1][3] =          0.50 * 0.75 ;

          trans[2][0] =                    0 ;
          trans[2][1] =                    0 ;
          trans[2][2] =                    1 ;
          trans[2][3] =                    0 ;

          trans[3][0] =                    0 ;
          trans[3][1] =                    0 ;
          trans[3][2] =                    0 ;
          trans[3][3] =                    1 ;

          for (ii=0; ii<4; ii++){

               for(jj=0; jj<4; jj++){

                    for(kk=0; kk<4; kk++){

                         trans[ii][jj] += trans[ii][kk] * prevtrans[kk][jj] ;

                    }
               }
          }

          prevtrans[0][0] =     trans[0][0] ;
          prevtrans[0][1] =     trans[0][1] ;
          prevtrans[0][2] =     trans[0][2] ;
          prevtrans[0][3] =     trans[0][3] ;
          prevtrans[1][0] =     trans[1][0] ;
          prevtrans[1][1] =     trans[1][1] ;
          prevtrans[1][2] =     trans[1][2] ;
          prevtrans[1][3] =     trans[1][3] ;
          prevtrans[2][0] =     trans[2][0] ;
          prevtrans[2][1] =     trans[2][1] ;
          prevtrans[2][2] =     trans[2][2] ;
          prevtrans[2][3] =     trans[2][3] ;
          prevtrans[3][0] =     trans[3][0] ;
          prevtrans[3][1] =     trans[3][1] ;
          prevtrans[3][2] =     trans[3][2] ;
          prevtrans[3][3] =     trans[3][3] ;

}

  printf("To close this program type 'quit' and hit the return key\n");
  printf("               \n");
  scanf("%d", &quit);

  return 0;

}


这是上面的C代码返回的最终矩阵:

0.4821  3.5870  11.68  381.22
0       1       0      76.875
0       0       5      0
0       0       0      5

最佳答案

这条线

                     trans[ii][jj] += trans[ii][kk] * prevtrans[kk][jj] ;


是不正确的。您仍在使用trans来计算结果矩阵时对其进行修改。您需要另一个矩阵来临时存储乘法结果。然后使用:

      // Store the resultant matrix in temp.
      for (ii=0; ii<4; ii++){

           for(jj=0; jj<4; jj++){

                temp[ii][jj] = 0.0;

                for(kk=0; kk<4; kk++){

                     temp[ii][jj] += trans[ii][kk] * prevtrans[kk][jj] ;

                }
           }
      }

      // Transfer the data from temp to trans
      for (ii=0; ii<4; ii++){

           for(jj=0; jj<4; jj++){
                trans[ii][jj] = temp[ii][jj];
           }
      }

关于c - 在C中将矩阵相乘,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/27844503/

相关文章:

C - 程序结构(避免全局变量、包含等)

r - 具有透明背景但字体可见的ggrepel标签

javascript - 旋转后平移缓冲区

c - C语言获取int后从用户处获取字符串

使用 getaddrinfo 连接套接字超时

python - R 和 Python 的输出值不同?

arrays - 如何将数组缩放到另一个长度以保存它在 R 中的近似值

python - 计算大型相关矩阵的内存有效方法?

c++ - 无法使用 lapack 构建 C++ 代码

c - 如何在 linux 中通过 exec 命令使用 wget 命令?