mzd_poly.h

Go to the documentation of this file.

#ifndef M4RIE_MZD_POLY_H
#define M4RIE_MZD_POLY_H

/******************************************************************************
*
*            M4RIE: Linear Algebra over GF(2^e)
*
*    Copyright (C) 2011 Martin Albrecht <martinralbrecht@googlemail.com>
*
*  Distributed under the terms of the GNU General Public License (GEL)
*  version 2 or higher.
*
*    This code is distributed in the hope that it will be useful,
*    but WITHOUT ANY WARRANTY; without even the implied warranty of
*    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
*    General Public License for more details.
*
*  The full text of the GPL is available at:
*
*                  http://www.gnu.org/licenses/
******************************************************************************/

#include <m4ri/m4ri.h>
#include "mzd_ptr.h"
#include "gf2x.h"
#include "blm.h"

typedef struct {
  mzd_t **x;   
  rci_t nrows; 
  rci_t ncols; 
  deg_t depth;   
} mzd_poly_t;

static inline mzd_poly_t *_mzd_poly_add(mzd_poly_t *C, const mzd_poly_t *A, const mzd_poly_t *B, unsigned int offset) {
  _mzd_ptr_add(C->x+offset, (const mzd_t**)A->x, (const mzd_t**)B->x, A->depth);
  return C;
}

static inline mzd_poly_t *mzd_poly_add(mzd_poly_t *C, const mzd_poly_t *A, const mzd_poly_t *B) {
  assert(C->depth >= A->depth && A->depth == B->depth);
  return _mzd_poly_add(C, A, B, 0);
}

static inline mzd_poly_t *mzd_poly_init(const deg_t d, const rci_t m, const rci_t n) {
  mzd_poly_t *A = (mzd_poly_t*)m4ri_mm_malloc(sizeof(mzd_poly_t));
  A->x = (mzd_t**)m4ri_mm_malloc(sizeof(mzd_t*)*(d+1));

  A->nrows = m;
  A->ncols = n;
  A->depth = d+1;

  for(int i=0; i<A->depth; i++)
    A->x[i] = mzd_init(m,n);
  return A;
}

static inline void mzd_poly_free(mzd_poly_t *A) {
  for(int i=0; i<A->depth; i++)
   mzd_free(A->x[i]);
  m4ri_mm_free(A->x);
  m4ri_mm_free(A);
}

static inline mzd_poly_t *_mzd_poly_adapt_depth(mzd_poly_t *A, const deg_t new_depth) {
  if (new_depth < A->depth) {
    for(int i=new_depth; i<A->depth; i++) {
      mzd_free(A->x[i]);
      A->x[i] = NULL;
    }
  } else {
    for(int i=A->depth; i<new_depth; i++) {
      A->x[i] = mzd_init(A->nrows,A->ncols);
    }
  }
  A->depth = new_depth;
  return A;
}

static inline mzd_poly_t *_mzd_poly_addmul_naive(mzd_poly_t *C, const mzd_poly_t *A, const mzd_poly_t *B) {
  if (C == NULL)
    C = mzd_poly_init(A->depth+B->depth-1, A->nrows, B->ncols);

  for(unsigned int i=0; i<A->depth; i++) {
    for(unsigned int j=0; j<B->depth; j++) {
      mzd_addmul(C->x[i+j], A->x[i], B->x[j], 0);
    }
  }
  return C;
}

static inline mzd_poly_t *_mzd_poly_addmul_karatsubs_balanced(mzd_poly_t *C, const mzd_poly_t *A, const mzd_poly_t *B) {
  assert(A->depth == B->depth);

  if (C == NULL)
    C = mzd_poly_init(A->depth+B->depth-1, A->nrows, B->ncols);
  switch(A->depth) {
  case 0:
    m4ri_die("depth 0: seriously?");
  case  1: mzd_addmul(C->x[0], A->x[0], B->x[0], 0); break;
  case  2:  _mzd_ptr_addmul_karatsuba2(NULL, C->x, (const mzd_t**)A->x, (const mzd_t**)B->x); break;
  case  3:  _mzd_ptr_addmul_karatsuba3(NULL, C->x, (const mzd_t**)A->x, (const mzd_t**)B->x); break;
  case  4:  _mzd_ptr_addmul_karatsuba4(NULL, C->x, (const mzd_t**)A->x, (const mzd_t**)B->x); break;
  case  5:  _mzd_ptr_addmul_karatsuba5(NULL, C->x, (const mzd_t**)A->x, (const mzd_t**)B->x); break;
  case  6:  _mzd_ptr_addmul_karatsuba6(NULL, C->x, (const mzd_t**)A->x, (const mzd_t**)B->x); break;
  case  7:  _mzd_ptr_addmul_karatsuba7(NULL, C->x, (const mzd_t**)A->x, (const mzd_t**)B->x); break;
  case  8:  _mzd_ptr_addmul_karatsuba8(NULL, C->x, (const mzd_t**)A->x, (const mzd_t**)B->x); break;
  case  9:  _mzd_ptr_addmul_karatsuba9(NULL, C->x, (const mzd_t**)A->x, (const mzd_t**)B->x); break;
  case 10: _mzd_ptr_addmul_karatsuba10(NULL, C->x, (const mzd_t**)A->x, (const mzd_t**)B->x); break;
  case 11: _mzd_ptr_addmul_karatsuba11(NULL, C->x, (const mzd_t**)A->x, (const mzd_t**)B->x); break;
  case 12: _mzd_ptr_addmul_karatsuba12(NULL, C->x, (const mzd_t**)A->x, (const mzd_t**)B->x); break;
  case 13: _mzd_ptr_addmul_karatsuba13(NULL, C->x, (const mzd_t**)A->x, (const mzd_t**)B->x); break;
  case 14: _mzd_ptr_addmul_karatsuba14(NULL, C->x, (const mzd_t**)A->x, (const mzd_t**)B->x); break;
  case 15: _mzd_ptr_addmul_karatsuba15(NULL, C->x, (const mzd_t**)A->x, (const mzd_t**)B->x); break;
  case 16: _mzd_ptr_addmul_karatsuba16(NULL, C->x, (const mzd_t**)A->x, (const mzd_t**)B->x); break;
  default:
    m4ri_die("Not implemented\n");
  }
  return C;
}

static inline mzd_poly_t *_mzd_poly_addmul_blm(mzd_poly_t *C, mzd_poly_t *A, mzd_poly_t *B, const blm_t *f) {
  assert(f!=NULL);
  assert(f->F->ncols == A->depth && f->G->ncols == B->depth);

  if (C == NULL)
    C = mzd_poly_init(A->depth+B->depth-1, A->nrows, B->ncols);

  _mzd_ptr_apply_blm(C->x, (const mzd_t**)A->x, (const mzd_t**)B->x, f);
  return C;
}


static inline mzd_poly_t *_mzd_poly_addmul_crt(mzd_poly_t *C, mzd_poly_t *A, mzd_poly_t *B) {
  int *p = crt_init(A->depth, B->depth);
  blm_t *f = blm_init_crt(NULL, A->depth, B->depth, p, 1);
  _mzd_poly_addmul_blm(C, A, B, f);
  blm_free(f);
  m4ri_mm_free(p);
  return C;
}

mzd_poly_t *_mzd_poly_addmul_ext1(mzd_poly_t *C, mzd_poly_t *A, mzd_poly_t *B);

static inline int mzd_poly_cmp(mzd_poly_t *A, mzd_poly_t *B) {
  int r = 0;
  if ((A->depth != B->depth) ) {
    if (A->depth < B->depth)
      return -1;
    else 
      return 1;
  }
  for(int i=0; i<A->depth; i++)
    r |= mzd_cmp(A->x[i],B->x[i]);
  return r;
}

static inline void mzd_poly_randomize(mzd_poly_t *A) {
  for(int i=0; i<A->depth; i++)
    mzd_randomize(A->x[i]);
}

#endif //M4RIE_MZD_POLY_H