src/xbt/xbt_matrix.c

   1 /* Copyright (c) 2006-2014. The SimGrid Team.
   2  * All rights reserved.                                                     */
   3
   4 /* This program is free software; you can redistribute it and/or modify it
   5  * under the terms of the license (GNU LGPL) which comes with this package. */
   6
   7 #include <stdio.h>
   8 #include "xbt/sysdep.h"
   9 #include "xbt/log.h"
  10 #include "xbt/matrix.h"
  11
  12 XBT_LOG_NEW_DEFAULT_SUBCATEGORY(xbt_matrix, xbt, "2D data storage");
  13
  14 /** @brief constructor */
  15 xbt_matrix_t xbt_matrix_new(int lines, int rows, const unsigned long elmsize, void_f_pvoid_t const free_f)
  16 {
  17   xbt_matrix_t res = xbt_new(s_xbt_matrix_t, 1);
  18   res->lines = lines;
  19   res->rows = rows;
  20   res->elmsize = elmsize;
  21   res->free_f = free_f;
  22   res->data = xbt_malloc(elmsize * lines * rows);
  23   return res;
  24 }
  25
  26 /** @brief Creates a matrix being a submatrix of another one */
  27 xbt_matrix_t xbt_matrix_new_sub(xbt_matrix_t from, int lsize, int rsize, int lpos, int rpos,
  28                                 pvoid_f_pvoid_t const cpy_f)
  29 {
  30   xbt_matrix_t res = xbt_matrix_new(lsize, rsize, from->elmsize, from->free_f);
  31   xbt_matrix_copy_values(res, from, lsize, rsize, 0, 0, lpos, rpos, cpy_f);
  32   return res;
  33 }
  34
  35 /** @brief destructor */
  36 void xbt_matrix_free(xbt_matrix_t mat)
  37 {
  38   if (mat) {
  39     if (mat->free_f) {
  40       for (unsigned i = 0; i < (mat->lines * mat->rows); i++) {
  41         mat->free_f((void *) &(mat->data[i * mat->elmsize]));
  42       }
  43     }
  44     free(mat->data);
  45     free(mat);
  46   }
  47 }
  48
  49 /** @brief Display the content of a matrix (debugging purpose)
  50  *  @param coords: boolean indicating whether we should add the coords of each cell to the output*/
  51 void xbt_matrix_dump(xbt_matrix_t matrix, const char *name, int coords, void_f_pvoid_t display_fun)
  52 {
  53   fprintf(stderr, ">>> Matrix %s dump (%u x %u)\n", name, matrix->lines, matrix->rows);
  54   for (unsigned i = 0; i < matrix->lines; i++) {
  55     fprintf(stderr, "  ");
  56     for (unsigned j = 0; j < matrix->rows; j++) {
  57       if (coords)
  58         fprintf(stderr, " (%u,%u)=", i, j);
  59       else
  60         fprintf(stderr, " ");
  61       display_fun(xbt_matrix_get_ptr(matrix, i, j));
  62     }
  63     fprintf(stderr, "\n");
  64   }
  65   fprintf(stderr, "<<< end_of_matrix %s dump\n", name);
  66 }
  67
  68 /** @brief Copy the values from the matrix src into the matrix dst
  69  * @param dst: destination
  70  * @param src: source
  71  * @param lsize: number of lines to copy
  72  * @param rsize: number of rows to copy
  73  * @param lpos_dst: line offset on destination matrix
  74  * @param rpos_dst: row offset on destination matrix
  75  * @param lpos_src: line offset on destination matrix
  76  * @param rpos_src: row offset on destination matrix
  77  * @param cpy_f: the function to use to copy the elements over
  78  */
  79 void xbt_matrix_copy_values(xbt_matrix_t dst, xbt_matrix_t src, unsigned int lsize, unsigned int rsize,
  80                             unsigned int lpos_dst, unsigned int rpos_dst,unsigned int lpos_src, unsigned int rpos_src,
  81                             pvoid_f_pvoid_t const cpy_f)
  82 {
  83   XBT_DEBUG ("Copy a %ux%u submatrix from %ux%u(of %ux%u) to %ux%u (of %ux%u)",
  84        lsize, rsize, lpos_src, rpos_src, src->lines, src->rows, lpos_dst, rpos_dst, dst->lines, dst->rows);
  85
  86   /* everybody knows that issue is between the chair and the screen (particularly in my office) */
  87   xbt_assert(src->elmsize == dst->elmsize);
  88   /* don't check free_f since the user may play weird games with this */
  89
  90   xbt_assert(lpos_src + lsize <= src->lines);
  91   xbt_assert(rpos_src + rsize <= src->rows);
  92
  93   xbt_assert(lpos_dst + lsize <= dst->lines);
  94   xbt_assert(rpos_dst + rsize <= dst->rows);
  95
  96   /* Lets get serious here */
  97   for (unsigned i = 0; i < rsize; i++) {
  98     if (cpy_f) {
  99       for (unsigned j = 0; j < lsize; j++)
 100         xbt_matrix_get_as(dst, j + lpos_dst, i + rpos_dst, void *) =
 101             cpy_f(xbt_matrix_get_ptr(src, j + rpos_src, i + lpos_src));
 102     } else {
 103       memcpy(xbt_matrix_get_ptr(dst, lpos_dst, i + rpos_dst), xbt_matrix_get_ptr(src, lpos_src, i + rpos_src),
 104              dst->elmsize * lsize);
 105     }
 106   }
 107 }
 108
 109 /** @brief Creates a new matrix of double filled with zeros */
 110 xbt_matrix_t xbt_matrix_double_new_zeros(int lines, int rows)
 111 {
 112   xbt_matrix_t res = xbt_matrix_new(lines, rows, sizeof(double), NULL);
 113
 114   memset(res->data, 0, res->elmsize * res->lines * res->rows);
 115   return res;
 116 }
 117
 118 /** @brief Creates a new matrix of double being the identity matrix */
 119 xbt_matrix_t xbt_matrix_double_new_id(int lines, int rows)
 120 {
 121   xbt_matrix_t res = xbt_matrix_double_new_zeros(lines, rows);
 122
 123   for (int i = 0; i < lines; i++)
 124     xbt_matrix_get_as(res, i, i, double) = 1;
 125   return res;
 126 }
 127
 128 /** @brief Creates a new matrix of double containing the sequence of numbers in order */
 129 xbt_matrix_t xbt_matrix_double_new_seq(int lines, int rows)
 130 {
 131   xbt_matrix_t res = xbt_matrix_new(lines, rows, sizeof(double), NULL);
 132
 133   for (int i = 0; i < lines * rows; i++)
 134     *(double *) &res->data[i * res->elmsize] = i;
 135
 136   return res;
 137 }
 138
 139 /** @brief add to C the result of A*B */
 140 void xbt_matrix_double_addmult(xbt_matrix_t A, xbt_matrix_t B, /*OUT*/ xbt_matrix_t C)
 141 {
 142   xbt_assert(A->lines == C->lines, "A->lines != C->lines (%u vs %u)", A->lines, C->lines);
 143   xbt_assert(B->rows == C->rows);
 144
 145   for (unsigned i = 0; i < C->lines; i++)
 146     for (unsigned j = 0; j < C->rows; j++)
 147       for (unsigned k = 0; k < B->lines; k++)
 148         xbt_matrix_get_as(C, i, j, double) += xbt_matrix_get_as(A, i, k, double) * xbt_matrix_get_as(B, k, j, double);
 149 }