/*
  Author:  Pate Williams (c) 1997

  Algorithm 1.4.3 (Order of an Element). See "A
  Course in Computational Agebraic Number Theory"
  by Henri Cohen page 25. Given a finite group
  G of cardinality h = |G|, and an element g in
  the group, this algorithm computes the order
  of g in G. This code uses trial division to
  factor the cardinality h. As an example
  consider the multiplicative group Z*19, where
  p in the program is 19. e is the order of g.
  -----
  g  e
  -----
  1  1
  18 2
  7  3
  8  6
  4  9
  2  18
  -----
*/

#include <ctype.h>
#include <stdio.h>
#include <stdlib.h>
#include "lip.h"

#define BOUND 10000000l

struct factor {int expon; verylong prime;};

void trial_division(long B, verylong zN, struct factor *f,
                    long *count)
{
  int e, one = 0, pri = 0;
  long p;
  verylong za = 0, zb = 0;

  *count = 0;
  zcopy(zN, &za);
  zpstart2();
  do {
    p = zpnext();
    if (zsmod(za, p) == 0) {
      e = 0;
      do {
        zsdiv(za, p, &zb);
        zcopy(zb, &za);
        e++;
      } while (zsmod(za, p) == 0);
      f[*count].expon = e;
      zintoz(p, &f[*count].prime);
      *count = *count + 1;
      one = zscompare(za, 1l) == 0;
      if (!one) pri = zprobprime(za, 5l);
    }
  } while (!one && !pri && p <= B);
  if (pri) {
    f[*count].expon = 1;
    zcopy(za, &f[*count].prime);
    *count = *count + 1;
  }
  if (p > B) {
    fprintf(stderr, "fatal error\ncan't factor number");
    exit(1);
  }
  zfree(&za);
  zfree(&zb);
}

void order(verylong zg, verylong zh,
           verylong zp, verylong *ze)
{
  long i, k;
  struct factor f[32];
  verylong za = 0, zq = 0, zr = 0, zg1 = 0;

  for (i = 0; i < 32; i++) f[i].prime = 0;
  trial_division(BOUND, zh, f, &k);
  zcopy(zh, ze);
  for (i = 0; i < k; i++) {
    zsexp(f[i].prime, f[i].expon, &za);
    zdiv(*ze, za, &zq, &zr);
    zcopy(zq, ze);
    zexpmod(zg, *ze, zp, &zg1);
    zcopy(f[i].prime, &za);
    while (zscompare(zg1, 1l) != 0) {
      zexpmod(zg1, za, zp, &zq);
      zcopy(zq, &zg1);
      zmul(*ze, za, &zq);
      zcopy(zq, ze);
    }
  }
  zfree(&za);
  zfree(&zq);
  zfree(&zr);
  zfree(&zg1);
  for (i = 0; i < 32; i++) zfree(&f[i].prime);
}

int main(void)
{
  char answer[256];
  verylong ze = 0, zg = 0, zh = 0, zp = 0;

  do {
    printf("g = "); zread(&zg);
    printf("p = "); zread(&zp);
    zsadd(zp, - 1l, &zh);
    order(zg, zh, zp, &ze);
    printf("e = "); zwriteln(ze);
    printf("another group (n or y)? ");
    scanf("%s", answer);
  } while (tolower(answer[0]) == 'y');
  zfree(&ze);
  zfree(&zg);
  zfree(&zh);
  zfree(&zp);
  return 0;
}