/*  This file constructs a PC presentation for a generalised
Nottingham group.
*/

#include "global.h"

#include <stdio.h>

extern int prime, rank, total;

/* Imports from group_elt.c */

void gcopy(group_elt f, group_elt g);
void initgroupelt(group_elt g);
void gprod(group_elt f, group_elt g, group_elt h);
void gpower(group_elt g, int n);
void groupinv(group_elt g, group_elt h);
int groupnotid(group_elt h);
void commutator(group_elt f, group_elt g, group_elt h);
void writegpelt(group_elt g);

/* Prototypes */

void normalform(group_elt g, vector v);
void printvec(vector v);
void powers(void);
void commutators(void);
void powers2(void);

void normalform(group_elt g, vector v)
{
   int i, j, k = 0 , m;
   group_elt h, b, c, inv;
   gcopy(h, g);
   for (i = rank + 1; i < total; i++)
   {
     for (j = 0; j < rank; j++)
     if (h[j][i])
     {
       v[k++] = h[j][i];
       initgroupelt(b);
       for (m = 0; m < rank; m++)
         b[m][rank - m] = 1; 
       b[j][i] = 1;
       groupinv(b, inv);
       gpower(inv , h[j][i]);
       gprod(inv, h, c);
       gcopy(h, c);
     }
     else v[k++] = 0;
   }
   if (k = groupnotid(h))
   {
     printf("Normal form failed, k = %d\n", k);
     exit(0);
   }
}

void printvec(vector v)
{
  int i, id = 1;
  for(i = 0; i < (total - rank -1) * rank; i++)
    if (v[i])
      {printf("a_%d^%d ", i + 1,  v[i]); id = 0;}
  if (id) printf(" 1");
  printf("\n");
}

void powers(void)
{
  int i, j, k, t = 1;
  vector v;
  group_elt g;
  for(i = rank + 1; i < total; i++)
    for(j = 0; j < rank; j++)   
    {
      initgroupelt(g);
      for(k = 0; k < rank; k++)
        g[k][rank - k] = 1;
      g[j][i] = 1;
      gpower(g, prime);
      normalform(g, v);
      printf("a_%d^%d = ",t++ , prime);
      printvec(v);
    }
}

void commutators(void)
{
  group_elt g, h, c;
  vector vec;
  int i, j, u, v, k, t, s = 0;
  for(i = rank + 1; i < total; i++)
    for(j = 0; j < rank; j++)
    {
      initgroupelt(g);
      for(k = 0; k < rank; k++)
        g[k][rank - k] = 1;
      g[j][i] = 1;
      s++;
      t = 0;
      for(u = rank + 1; u < total; u++)
        for(v = 0; v < rank; v++)
        {
         t++;
         if (t>s)
         {
           initgroupelt(h);
           for(k = 0; k < rank; k++)
             h[k][rank - k] = 1;
           h[v][u] = 1;
           commutator(h, g, c);
           normalform(c, vec); 
           printf("[a_%d,a_%d]=",t,s);          
           printvec(vec);
           printf("\n");
         }
        }
    }
}       

void powers2(void)
{
  int i, j, k, t = 1;
  vector *pv;
  group_elt g;
  char c[2];
  for(i = rank + 1; i < total; i++)
    for(j = 0; j < rank; j++)   
    {
      initgroupelt(g);
      for(k = 0; k < rank; k++)
        g[k][rank - k] = 1;
      g[j][i] = 1;
      writegpelt(g);
      gpower(g, prime);
      printf("p-th power is \n");
      writegpelt(g);
      printf("more? 'y' or 'n'\n");
      scanf("%s", c);
      if (c[0] != 'y') exit(0);
    }
}