// gcc -O3 -o pi pi.c
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <sys/mman.h>
#include <sys/types.h>
#include <sys/wait.h>
#include <sys/time.h> 
#include <unistd.h>

#define MAX_HEXDIG 15  // maximum number of hexadecimal digits allowed
#define MAX_PR 8       // maximum number or processes allowed

int digits=2;
int procs=1;
long long unsigned int t;    // fator de escala para o valor 1
long long unsigned int iter; // numero de iteracoes

// calcula pi segundo Leibniz usando n digitos fracionarios hexadecimais"
long long unsigned int pi(long long unsigned int a) { // argumento: valor inicial do denomidar na sequencia
  long long unsigned int m,     // valor 4, ja com fator de escala
                         d,     // difference to next value of pi
                         p=0;   // o resultado, em escala de n digitos hexadecimais
  int b; // incremento para a
  m=4*t;
  b=4*procs;
                         
  while(1) { 
    d = m/a - m/(a+2); // dois elementos da sequencia de Leibniz
    if(d==0) break;  // termine quando pi nao varia mais de uma rodada para a proxima
    p = p + d;
    a = a+b;
  }
  iter = a/b;
  return p;
}

static struct { long long unsigned int pi, iter; } *pr; // parts of pi and iter

int main(int argc, char *argv[]) {
  static long long int p; // resultado
  struct timeval start_time, end_time;
  double run_time;
  int i;  // para os lacos for
  char p_exact[] = "3.243F6A8885A308D313";
  char  p_str[MAX_HEXDIG+3]; // resultado na forma de string

  if(argc>1) digits = atoi(argv[1]);
  if(digits<1 || digits>MAX_HEXDIG) {
    fprintf(stderr,"Can use at most %d digits, not %d\n.", MAX_HEXDIG, digits);
    exit(1);
  }
  t=(long long unsigned int)1<<(4*digits);
  if(argc>2) procs = atoi(argv[2]);
  if(procs<1 || procs>MAX_PR) {
    fprintf(stderr,"Can use launch most %d processes, not %d\n.", MAX_PR, procs);
    exit(1);
  }

  pr = mmap(NULL, (MAX_PR-1)*(sizeof *pr), PROT_READ | PROT_WRITE, 
            MAP_SHARED | MAP_ANONYMOUS, -1, 0);

  gettimeofday(&start_time, NULL);

  for(i=0; i<procs-1; i++)
    if (fork() == 0) {     // run another seperate process
      pr[i].pi = pi(1+4*(i+1));
      pr[i].iter = iter;
      exit(EXIT_SUCCESS);
    }
 
  p = pi(1);        // run first process
  for(i=0; i<procs-1; i++) wait(NULL); // wait for the others to finish
  for(i=0; i<procs-1; i++) {
    p = p + pr[i].pi;
    iter += pr[i].iter;
  }  

  gettimeofday(&end_time, NULL);
  run_time = (double) (end_time.tv_usec - start_time.tv_usec) / 1000000 +
             (double) (end_time.tv_sec  - start_time.tv_sec);

  munmap(pr, (MAX_PR-1)*(sizeof *pr));

  sprintf(p_str+1, "%lX", p); // put fracional separator (point) after first digit
  p_str[0] = p_str[1];
  p_str[1] = '.';
  p_exact[2+digits] = '\0';
  printf("resultado   em hexadecimal: %s\n", p_str);
  printf("valor exato em hexadecimal: %s\n", p_exact);
  printf("resultado   em decimal: %1.15f\n", ((double)p)/t);
  printf("valor exato em decimal: 3.141592653589793\n");
  printf("tempo de execucao: %.1f segundos\n", run_time);
  printf("%ld iteracoes, %d ns/iteracao", iter, (int)(run_time*1000000000/(iter/procs)));
  if(procs>1) printf(", %d iteracoes em parallelo\n", procs);
  else        printf("\n");
  return 0;
}