c/c++语言开发共享C语言实现常见的矩阵运算函数

1.矩阵转置函数

  void matrix_t(double **a_matrix, const double **b_matrix, int krow, int kline)  ////////////////////////////////////////////////////////////////////////////  //  a_matrix:转置后的矩阵  //  b_matrix:转置前的矩阵  //  krow    :行数  //  kline   :列数  ////////////////////////////////////////////////////////////////////////////  {      int k, k2;           for (k = 0; k < krow; k++)      {          for(k2 = 0; k2 < kline; k2++)          {              a_matrix[k2][k] = b_matrix[k][k2];          }      }  }  

2.矩阵加（减）法函数

  void matrix_a(double **a_matrix, const double **b_matrix, const double **c_matrix,                       int krow, int kline, int ktrl)  ////////////////////////////////////////////////////////////////////////////  //  a_matrix=b_matrix+c_matrix  //   krow   :行数  //   kline  :列数  //   ktrl   :大于0: 加法  不大于0:减法  ////////////////////////////////////////////////////////////////////////////  {      int k, k2;        for (k = 0; k < krow; k++)      {          for(k2 = 0; k2 < kline; k2++)          {              a_matrix[k][k2] = b_matrix[k][k2]                  + ((ktrl > 0) ? c_matrix[k][k2] : -c_matrix[k][k2]);           }      }  }

3.矩阵乘法函数

  void matrix_m(double **a_matrix, const double **b_matrix, const double **c_matrix,                  int krow, int kline, int kmiddle, int ktrl)  ////////////////////////////////////////////////////////////////////////////  //  a_matrix=b_matrix*c_matrix  //  krow  :行数  //  kline :列数  //  ktrl  : 大于0:两个正数矩阵相乘 不大于0:正数矩阵乘以负数矩阵  ////////////////////////////////////////////////////////////////////////////  {      int k, k2, k4;      double stmp;        for (k = 0; k < krow; k++)           {          for (k2 = 0; k2 < kline; k2++)             {              stmp = 0.0;              for (k4 = 0; k4 < kmiddle; k4++)                {                  stmp += b_matrix[k][k4] * c_matrix[k4][k2];              }              a_matrix[k][k2] = stmp;          }      }      if (ktrl <= 0)         {          for (k = 0; k < krow; k++)          {              for (k2 = 0; k2 < kline; k2++)              {                  a_matrix[k][k2] = -a_matrix[k][k2];              }          }      }  }

4.矩阵求逆函数

  int  matrix_inv(double **a_matrix, int ndimen)  ////////////////////////////////////////////////////////////////////////////  //  a_matrix:矩阵  //  ndimen :维数  ////////////////////////////////////////////////////////////////////////////  {      double tmp, tmp2, b_tmp[20], c_tmp[20];      int k, k1, k2, k3, j, i, j2, i2, kme[20], kmf[20];      i2 = j2 = 0;        for (k = 0; k < ndimen; k++)        {          tmp2 = 0.0;          for (i = k; i < ndimen; i++)            {              for (j = k; j < ndimen; j++)                {                  if (fabs(a_matrix[i][j] ) <= fabs(tmp2))                       continue;                  tmp2 = a_matrix[i][j];                  i2 = i;                  j2 = j;              }            }          if (i2 != k)           {              for (j = 0; j < ndimen; j++)                 {                  tmp = a_matrix[i2][j];                  a_matrix[i2][j] = a_matrix[k][j];                  a_matrix[k][j] = tmp;              }          }          if (j2 != k)           {              for (i = 0; i < ndimen; i++)                {                  tmp = a_matrix[i][j2];                  a_matrix[i][j2] = a_matrix[i][k];                  a_matrix[i][k] = tmp;              }              }          kme[k] = i2;          kmf[k] = j2;          for (j = 0; j < ndimen; j++)            {              if (j == k)                 {                  b_tmp[j] = 1.0 / tmp2;                  c_tmp[j] = 1.0;              }              else               {                  b_tmp[j] = -a_matrix[k][j] / tmp2;                  c_tmp[j] = a_matrix[j][k];              }              a_matrix[k][j] = 0.0;              a_matrix[j][k] = 0.0;          }          for (i = 0; i < ndimen; i++)            {              for (j = 0; j < ndimen; j++)                {                  a_matrix[i][j] = a_matrix[i][j] + c_tmp[i] * b_tmp[j];              }            }      }      for (k3 = 0; k3 < ndimen;  k3++)         {          k  = ndimen - k3 - 1;          k1 = kme[k];          k2 = kmf[k];          if (k1 != k)             {              for (i = 0; i < ndimen; i++)                {                  tmp = a_matrix[i][k1];                  a_matrix[i][k1] = a_matrix[i][k];                  a_matrix[i][k] = tmp;              }            }          if (k2 != k)             {              for(j = 0; j < ndimen; j++)                {                  tmp = a_matrix[k2][j];                  a_matrix[k2][j] = a_matrix[k][j];                  a_matrix[k][j] = tmp;              }          }      }      return (0);  }  

5.矩阵乔里斯基分解函数

  void chol(double **a_matrix, const double **b_matrix, int ndimen)  ////////////////////////////////////////////////////////////////////////////  //  输入参数:  //      b_matrix:  对称正定方阵    ndimen: 矩阵维数  //  返回值:  //      a_matrix: 下三角矩阵  ////////////////////////////////////////////////////////////////////////////  {      int i, j, r;      double m = 0;         static double **c_matrix;      static int flag = 0;        if (flag == 0)      {          flag = 1;          c_matrix = (double **)malloc(ndimen * sizeof(double *));            for (i = 0; i < ndimen; i++)              c_matrix[i] = (double *)malloc(ndimen * sizeof(double));      }        for (i = 0; i < ndimen; i++)      {          for (j = 0; j < ndimen; j++)               c_matrix[i][j] = 0;      }        c_matrix[0][0] = sqrt(b_matrix[0][0]);        for (i = 1; i < ndimen; i++)      {          if (c_matrix[0][0] != 0)               c_matrix[i][0] = b_matrix[i][0] / c_matrix[0][0];      }        for (i = 1; i < ndimen; i++)      {          for (r = 0; r < i; r++)      m = m + c_matrix[i][r] * c_matrix[i][r];            c_matrix[i][i] = sqrt(b_matrix[i][i] - m);          m = 0.0;            for (j = i + 1; j < ndimen; j++)          {              for (r = 0; r < i; r++)      m = m + c_matrix[i][r] * c_matrix[j][r];              c_matrix[j][i] = (b_matrix[i][j] - m) / c_matrix[i][i];              m = 0;          }      }        for (i = 0; i < ndimen; i++)      {          for (j = 0; j < ndimen; j++)               a_matrix[i][j] = c_matrix[i][j];      }  }