Ramón Casero Cañas on Mon, 20 Jan 2003 10:49:56 +0100 (MET)


[Date Prev] [Date Next] [Thread Prev] [Thread Next] [Date Index] [Thread Index]

Re: how to free pari's memory


Karim BELABAS wrote:

> On Sat, 18 Jan 2003, Ramón Casero Cañas wrote:
>
>> Karim BELABAS wrote:
>>
>>
>>> On Thu, 16 Jan 2003, Ramón Casero Cañas wrote:
>>
>>
>>> Anyway, if you want to use PARI, then remove it completely [ not a good
idea
>>> if you do this more than once : it involves a very significant overhead
],
>>> you have to do it by hand:
>>>
>>> void freeall(void);
>>> void pari_sig_init(void (*f)(int));
>>>
>>> freeall();
>>> if (INIT_SIG) pari_sig_init(SIG_DFL);


I get a

g++ -o lib/gamatrix.o -c -Wall   -ftemplate-depth-50 gamatrix.cpp
gamatrix.cpp: In function `bool inverse_mat_mod2(bool **, unsigned int,
bool)':
gamatrix.cpp:82: `INIT_SIG' undeclared (first use this function)
gamatrix.cpp:82: (Each undeclared identifier is reported only once
gamatrix.cpp:82: for each function it appears in.)
make: *** [lib/gamatrix.o] Error 1

at compilation. Code's:

#include <pari/pari.h>
void freeall(void); // thanks, Karim
void pari_sig_init(void (*f)(int));  // thanks, Karim

/* inverse_mat_mod2: overwrites a matrix modulo 2 with its inverse;
 * input is a bool matrix, but internally, PARI library types are
 * used */
bool inverse_mat_mod2(bool **_a, unsigned _size) {
  return inverse_mat_mod2(_a, _size, true);
}

bool inverse_mat_mod2(bool **_a, const unsigned _size, const bool _CHECK) {
  GEN a;           // temp matrix by columns
  unsigned i, j;   // counter
  unsigned sizepari = _size + 1;
  GEN z;           // auxiliary element of the matrix
  GEN v;           // auxiliary element of the matrix
  long paristacksize = 500000;
  long total_size = sizepari * sizepari + _size * _size * 9; // thanks,
Karim

  paristacksize = sizeof(long) * total_size * 10; // thanks again
  pari_init(paristacksize, 2); // init PARI stack

  /* allocate memory for the temp matrix
   * note: remember that PARI matrices of mxn need (m+1)x(n+1), and
   * that elements go a[1], a[2], ..., while a[0] is for type
   * information */
  std::cout << "inverse_mat_mod2: size = " << _size << std::endl;
  a = cgetg(sizepari, t_MAT); // `a' has _size columns
  for (i = 1; i <= _size; ++i) { // allocate `a' columns of _size elements
    a[i] = lgetg(sizepari, t_COL);
  }
  for (i = 1; i <= _size; ++i) { // allocate space for elements
    v = (GEN)a[i];
    for (j = 1; j <= _size; ++j) {
      v[j] = lgetg(3, t_INTMOD);
      z = (GEN)v[j];
      z[1] = lstoi(2); // modulo 2
      if (_a[j-1][i-1]) z[2] = un; else z[2] = zero;//lstoi(0);
    }
  }
  std::cout << "inverse_mat_mod2: assigment done" << std::endl;

  /* invert matrix */
  GEN ainv;
  ainv = ginv(a);

  /* check whether the inverse is valid */
  if (_CHECK) {
    GEN a_ainv;
    a_ainv = gmul(a, ainv);
    for (i = 1; i <= _size; ++i) {
      v = (GEN)a_ainv[i];
      for (j = 1; j <= _size; ++j) {
        z = (GEN)v[j];
        if ((i == j) != itos((GEN)z[2])) return false;
      }
    }
  }

  /* for no _CHECK or valid inversion, overwrite input array with
     inverse matrix */
  for (i = 1; i <= _size; ++i) {
    v = (GEN)ainv[i];
    for (j = 1; j <= _size; ++j) {
      z = (GEN)v[j];
      _a[j-1][i-1] = (bool)itos((GEN)z[2]);
    }
  }

  /* free PARI's memory (thanks Karim) */
  freeall();
  if (INIT_SIG) pari_sig_init(SIG_DFL);

  return true;
}


-- 
+++ GMX - Mail, Messaging & more  http://www.gmx.net +++
NEU: Mit GMX ins Internet. Rund um die Uhr für 1 ct/ Min. surfen!