/*
Copyright (C) 1998, 2000 Harald A. Helfgott

This program is free software; you can redistribute it and/or modify it 
under the terms of the GNU General Public License as published by the
Free Software Foundation; either version 2, or (at your option)
any later version.

This library 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.

If you do not have a copy of the GNU General Public License, write to the Free
Software Foundation, 59 Temple Place - Suite 330, Boston, MA  02111-1307, 
USA. */
#ifdef INF
#include <gmp.h>
#endif

#ifdef LDOUB
typedef long double num;
#else
#ifdef DOUB
typedef double num;
#else
typedef float num;
#endif
#endif

typedef struct {
#ifdef INF
  mpq_t q;
#else
  num q;
#endif
  int i;
} couple;

typedef struct {
  couple *a;    
  int memside;
  int i0,j0;
  int side;
} matr;

typedef struct {
  char *a;
  int memside;
  int i0,j0;
  int side;
} chatr;

#define el(mat,i,j) ((mat).a+(mat).j0+(j)+((mat).i0+(i))*((mat).memside))
#define chel(mat,i,j) ((mat).a+(mat).j0+(j)+((mat).i0+(i))*((mat).memside))
void crap_out();
#ifdef INF
void coup_set(couple *a, int q);
#else
void coup_set(couple *a, num q); /* sets epsilon (i) to zero, too */
#endif
void coup_e_set(couple *a, int i);
void coup_ass(couple *a, couple b);
void coup_add(couple *a, couple b, couple c);
void coup_sub(couple *a, couple b, couple c);
void coup_mul(couple *a, couple b, couple c);
void coup_div(couple *a, couple b, couple c);
/* all of these allow a and b, or a and c, to be the same */

void coup_init(couple *a); /*initializes couple.q, sets couple.i to zero */
void coup_kill(couple *a);
size_t coup_in(FILE *stream, int base, couple *op);
size_t coup_out(FILE *stream, int base, couple op);    

void matr_in(FILE *stream, int base, matr *mat);
void matr_out(FILE *stream, int base, matr mat, int offset);
void matr_create(matr *mat, int memside);
void matr_cpy(matr *mat, matr mat2);
void matr_kill(matr *mat);

void grow(matr *neu, matr weight, matr alt);
/* Pre: neu->memside, weight.memside, weight.side even */
/* neu is taken to be allocated and initialized */
/* Post: neu->side = weight.side = alt.side+2*/
void renew(matr *neu, matr alt, couple *fact);
/* Pre: neu->memside, alt.memside, alt.side even */
/* neu and fact are taken to be allocated and initialized */
/* Post: neu->side = alt.side-2 */

void forward(matr *neu, matr alt, couple *fact, int times);
/* renew alt times times */
/* beware: times=0 will lead to a copy of alt sharing the matrix with alt */
void backstep(matr *ar, int numdigs, int p);
matr edgprob(matr weight, couple *fact);
matr logedgprob(matr weight, couple *fact);