matrix/matrix.c

   1 /* matrix.c -- basic matrix operations
   2  *
   3  * author: hackbard@hackdaworld.dyndns.org
   4  *
   5  */
   6
   7 #include "matrix.h"
   8
   9 int matrix_init(t_matrix *matrix,int outfd) {
  10
  11  dprintf(outfd,"[matrix] initializing matrix api ...\n");
  12
  13  matrix->outfd=outfd;
  14  matrix->type=0;
  15  matrix->m=0;
  16  matrix->n=0;
  17  matrix->z=0;
  18  matrix->m1=NULL;
  19  matrix->m2=NULL;
  20
  21  return M_SUCCESS;
  22 }
  23
  24 int matrix_alloc_data(t_matrix *matrix) {
  25
  26   /* there is always one matrix */
  27   if((matrix->m==0)||(matrix->n==0)) {
  28     dprintf(matrix->outfd,"[matrix] wrong dimensions ...\n");
  29     return M_DIM_FAILURE;
  30   }
  31
  32   matrix->m1=malloc(matrix->m*matrix->n*sizeof(double));
  33   if(matrix->m1==NULL) {
  34     dprintf(matrix->outfd,"[matrix] malloc failed\n");
  35     return M_ALLOC_FAIL;
  36   }
  37
  38   /* but often there are more ... */
  39   if(matrix->type&MULT) {
  40     if(matrix->z==0) {
  41       dprintf(matrix->outfd,"[matrix] wrong dimensions ...\n");
  42       return M_DIM_FAILURE;
  43     }
  44
  45     matrix->m2=malloc(matrix->n*matrix->z*sizeof(double));
  46     if(matrix->m2==NULL) {
  47       dprintf(matrix->outfd,"[matrix] malloc failed\n");
  48       return M_ALLOC_FAIL;
  49     }
  50     matrix->res=malloc(matrix->m*matrix->z*sizeof(double));
  51     if(matrix->res==NULL) {
  52       dprintf(matrix->outfd,"[matrix] malloc failed\n");
  53       return M_ALLOC_FAIL;
  54     }
  55   } else {
  56     dprintf(matrix->outfd,"[matrix] type of operation not supported by now\n");
  57     return M_NOT_SUPPORTED;
  58   }
  59
  60   return M_SUCCESS;
  61 }
  62
  63 int matrix_shutdown(t_matrix *matrix) {
  64
  65   dprintf(matrix->outfd,"[matrix] shutdown\n");
  66
  67   if(matrix->m1!=NULL) free(matrix->m1);
  68   if(matrix->m2!=NULL) free(matrix->m2);
  69   if(matrix->res!=NULL) free(matrix->res);
  70
  71   return M_SUCCESS;
  72 }
  73
  74 /* matrix operations */
  75
  76 int matrix_calc(t_matrix *matrix) {
  77
  78   unsigned int i,j,k,off,off1,off2;
  79
  80   if(!(matrix->type&MULT)) {
  81     dprintf(matrix->outfd,"[matrix] only MULT supported by now ...\n");
  82     return M_NOT_SUPPORTED;
  83   }
  84
  85   /* mult operation, an O(N^3) problem
  86    *
  87    * maybe speedup by different order of mult?
  88    *
  89    */
  90   for(i=0;i<matrix->z;i++) {
  91     off1=0;
  92     off=i;
  93     for(j=0;j<matrix->n;j++) {
  94       off+=matrix->z;
  95       off2=i;
  96       matrix->res[off]=0;
  97       for(k=0;k<matrix->m;k++) {
  98         matrix->res[off]+=(matrix->m1[off1]*matrix->m2[off2]);
  99         off1+=1;
 100         off2+=matrix->z;
 101       }
 102     }
 103   }
 104
 105   return M_SUCCESS;
 106 }