#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include "pari.h"

#if defined (USG) || defined (__SVR4) || defined (_UNICOS) || defined(HPUX)
#include <time.h>

int
cputime ()
{
  if (CLOCKS_PER_SEC < 100000)
    return clock () * 1000 / CLOCKS_PER_SEC;
  return clock () / (CLOCKS_PER_SEC / 1000);
}
#else
#include <sys/types.h>
#include <sys/resource.h>

int
cputime ()
{
  struct rusage rus;

  getrusage (0, &rus);
  return rus.ru_utime.tv_sec * 1000 + rus.ru_utime.tv_usec / 1000;
}
#endif

GEN Fp_div(GEN,GEN,GEN);
GEN red_montgomery(GEN T, GEN N, ulong inv);

int
main(int argc, char *argv[])
{
  int n, prec, st, st2, N, i;
  GEN x, y, xmody, amax, z, r, s, X, mont;
  pari_sp ltop;
  ulong inv;

  if (argc != 2 && argc != 3) {
    fprintf(stderr, "Usage: timing digits [N]\n"); exit(1);
  }
  printf("%s\n",PARIVERSION);
  n = atoi(argv[1]);
  if (argc==3)  N = atoi(argv[2]);
  prec = (int) ( n * log(10.0) / log(2.0) + 1.0 );
  printf("prec=%u\n", n);

  pari_init(10000000, 2);

  x = randomi(int2n(prec));
  y = randomi(int2n(prec));
  X = randomi(int2n(2 * prec));
  amax = sqrti(shifti(y,-1));
  for(;;) {
    GEN a = randomi(amax), g = gcdii(a,y);
    if (equali1(g)) {
      xmody = Fp_div(randomi(amax), a, y);
      break;
    }
  }
  mont = x; if (!mod2(mont)) mont = addis(mont,1); /* odd */
  inv = (ulong) -invmod2BIL(mod2BIL(mont));

#define TIME(__s__, __EXPR__) \
  N=1;  st = cputime(); \
  do { \
    for (i=0;i<N;i++) __EXPR__, avma = ltop; \
    N=2*N; \
    st2=cputime(); \
  } while (st2-st<1000); \
  printf("%-12s took %f ms (%d eval in %d ms)\n", \
         __s__,(double)(st2-st)/(N-1),N-1,st2-st);

  ltop = avma;
  TIME("addii", z = addii (x, y));
  TIME("sqri", z = sqri (x));
  TIME("mulii", z = mulii (x, y));
  TIME("sqrti", z = sqrti (x));

  TIME("diviiexact", z = diviiexact (X, y));
  TIME("divii", z = divii (X, y));
  TIME("dvmdii", z = dvmdii (X, y, &r));
  TIME("redc", z = red_montgomery(X, mont, inv));

  TIME("gcdii", z = gcdii (x, y));
  TIME("bezout", z = bezout (x, y, &r,&s));
  TIME("invmod", (void)invmod (x, y, &s));
  TIME("ratlift", (void)Fp_ratlift(xmody, y, amax,amax,&r,&s));
  return 0;
}