#include <stdio.h>
#include <time.h>

#define OVERFLOW	1
#define CASAS		1100
#define RESUMO		1000
#define CESQUERDA	3
#define BASE		7
#define FORMATO		"%01d"

#define MAIS  		'+'
#define MENOS 		'-'


/* A estrutura do numero de ponto flutuante foi fixada em dois digitos
 * antes do ponto e CASAS digitos apos o ponto e um espaco para terminar
 * o string.
 */
typedef struct
	{
	char sinal;
	int digitos[CASAS + 1 + CESQUERDA];
	} numero;


int somaxy (numero *resultado, numero *x, numero *y)
	{
	numero *p;
	int soma;
	int acumulador;
	int i;

	if ((x->sinal) == (y->sinal))
		{
		resultado->sinal = x->sinal;

		acumulador = 0;
		for (i = CASAS + CESQUERDA; i >= 0; i--)
			{
			soma = x->digitos[i] + y->digitos[i] + acumulador;
			resultado->digitos[i] = soma % BASE;
			acumulador = soma / BASE;
			};
		}
	else
		{
		/* verifica qual numero e maior */
		for (i = 0; i <= CASAS + CESQUERDA; i++)
			if (x->digitos[i] != y->digitos[i])
				break;

		/* por elegancia faz o sinal do resultado ser +0 caso x e y
		 * sejam iguais.
		 */
		if (i > CASAS + CESQUERDA)
			{
			x->sinal = MAIS;
			y->sinal = MAIS;
			};

		/* se y for maior que x entao troca x com y */
		if (y->digitos[i] > x->digitos[i])
			{
			p = x;
			x = y;
			y = p;
			};

		/* agora que x >= y faz a subtracao */
		resultado->sinal = x->sinal;

		/* faz o truque de emprestar BASE a priori para facilitar o algoritmo.
		 * este truque e possivel porque nunca havera a necessidade de tomar
		 * mais que BASE emprestado.
		 */
		acumulador = 0;
		for (i = CASAS + CESQUERDA; i >= 0; i--)
			{
			soma = x->digitos[i] - y->digitos[i] + acumulador + BASE;
			resultado->digitos[i] = soma % BASE;
			acumulador = soma / BASE - 1;
			};
		};

	if (acumulador != 0)
		return (OVERFLOW);

	return (0);
	}


int mult[2 * CASAS + 2 * CESQUERDA + 2];

int multxy (numero *resultado, numero *x, numero *y)
	{
	long soma;
	long acumulador;
	int yi;
	int i, j, k;

	/* inicializa a variavel com o valor completo de x vezes y */
	for (i = 0; i <= 2 * CASAS + 2 * CESQUERDA + 1; i++)
		mult[i] = 0;

	for (i = CASAS + CESQUERDA; i >= 0; i--)
		{
		acumulador = 0;
		for (j = CASAS + CESQUERDA; j >= 0; j--)
			{
			k = i + j + 1 + CESQUERDA;

			soma = (long)mult[k] + (long)(x->digitos[i]) * (long)(y->digitos[j]) + acumulador;
			mult[k] = soma % BASE;
			acumulador = soma / BASE;
			};

		k--;
		while ((k >= 0) && (acumulador != 0))
			{
			soma = mult[k] + acumulador;
			mult[k] = soma % BASE;
			acumulador = soma / BASE;
			k--;
			};

		/* e impossivel dar overflow neste ponto mas o codigo esta ai para
		 * desencargo de conciencia
		 */
		if (acumulador > 0)
			return (OVERFLOW);
		};

	if (x->sinal == y->sinal)
		resultado->sinal = MAIS;
	else
		resultado->sinal = MENOS;

	for (i = 0; i <= CASAS + CESQUERDA; i++)
		resultado->digitos[i] = mult[i + 2 * CESQUERDA + 1];

	/* agora sim faz o teste de overflow pois o algoritmo pressupoe apenas
	 * CESQUERDA casa do lado esquerdo da virgula
	 */
	if (mult[0] != 0)
		return (OVERFLOW);

	return (0);
	}


int divxn (numero *resultado, numero *x, long n)
	{
	long dividendo;
	char sinal;
	int i;

	/* testa se o valor de n nao vai ocasionar overflow */
	if ((n == 0) || (abs (n) > 200000))
		return (OVERFLOW);

	if (n > 0)
		sinal = MAIS;
	else
		sinal = MENOS;

	if (x->sinal == sinal)
		resultado->sinal = MAIS;
	else
		resultado->sinal = MENOS;

	dividendo = 0;
	for (i = 0; i <= CASAS + CESQUERDA; i++)
		{
		dividendo = BASE * dividendo + x->digitos[i];

		resultado->digitos[i] = dividendo / n;
		dividendo = dividendo % n;
		};

	return (0);
	}


int ezero (numero *x)
   {
   int i;

   for (i = 0; i <= CASAS + CESQUERDA; i++)
      if (x->digitos[i] != 0)
	 return (1);

   x->sinal = MAIS;
   return (0);
   }


int printx (numero *x)
	{
	int i;

	putchar (x->sinal);
	for (i = 0; i <= CASAS + CESQUERDA; i++)
		{
		printf (FORMATO, x->digitos[i]);
		if (i == 0 + CESQUERDA)
			putchar ('.');
		};
	putchar ('\n');

	return (0);
	}


int loadx (numero *x)
	{
	FILE *src;
	int valor;
	int i, j;
	int ch;

	x->sinal = '+';
	for (i = 0; i <= CASAS + CESQUERDA; i++)
		x->digitos[i] = 0;
	x->digitos[0 + CESQUERDA] = 3;

	src = fopen ("pi.dat", "rt");
	if (src == NULL)
		return (0);

	ch = fgetc (src);
	if ((ch != '+') && (ch != '-'))
		{
		fclose (src);
		return (0);
		};

	x->sinal = ch;
	for (i = 0; i <= CASAS + CESQUERDA; i++)
		{
		valor = 0;

		for (j = BASE; j > 1; j = j / 10)
			{
			ch = fgetc (src);
			if ((ch < '0') || (ch > '9'))
				break;
			valor = valor * 10 + (ch - '0');
			};

		x->digitos[i] = valor;
		};

   fclose (src);
   return (0);
   }


int savex (numero *x)
	{
	FILE *dst;
	int i;

	dst = fopen ("pi.dat", "wt");
	if (dst == NULL)
		return (0);

	fputc (x->sinal, dst);
	for (i = 0; i <= CASAS + CESQUERDA; i++)
		fprintf (dst, FORMATO, x->digitos[i]);

	fclose (dst);
	return (0);
	}


int main ()
	{
	numero x, x2, fator, pi, tmp;
	long n;
	int i, j;
	float segundos;
	long contador[BASE];

	loadx (&pi);
	i = 0;
	n = 0;
	do
		{
		printf ("Pi(%d): ", ++i);
                printx (&pi);

		segundos = (float)clock() / (float)CLOCKS_PER_SEC;
		printf ("Iteracao %d: %d Casas Fracionarias na Base %d em %04.2f segundos (nmax=%ld)\n",
	    	i - 1, CASAS, BASE, segundos, n);

		x.sinal = pi.sinal;
		for (j = 0; j <= CASAS + CESQUERDA; j++)
			x.digitos[j] = pi.digitos[j];

		multxy (&x2, &x, &x);

		fator.sinal = x.sinal;
		for (j = 0; j <= CASAS + CESQUERDA; j++)
			fator.digitos[j] = x.digitos[j];

		n = 1;
		while (ezero (&fator) != 0)
			{
			multxy (&fator, &fator, &x2);

			n++;
			divxn (&fator, &fator, n);
			n++;
			divxn (&fator, &fator, n);

			if (fator.sinal == MAIS)
				fator.sinal = MENOS;
			else
				fator.sinal = MAIS;

			somaxy (&x, &x, &fator);
			};

		somaxy (&pi, &pi, &x);
		savex (&pi);
		}
	while ((ezero (&x) != 0) && (i < 10));

	printf ("\n");
	printf ("Pi(%d): ", i);
	printx (&pi);
	
	segundos = (float)clock() / (float)CLOCKS_PER_SEC;
	printf ("Iteracao %d: %d Casas Fracionarias na Base %d em %04.2f segundos\n",
   		i, CASAS, BASE, segundos);

	for (j = 0; j < BASE; j++)
	  contador[j] = 0;
	
	for (i = 0; i <= CASAS + CESQUERDA; i++)
	  for (j = 0; j < BASE; j++)
	    if (pi.digitos[i] == j)
	      contador[j]++;
	
	printf ("\n");
	printf ("Estatistica: (%d casas, com %d casas inteiras e %d casas fracionarias\n", CASAS + CESQUERDA, CESQUERDA, CASAS);
	
	for (j = 0; j < BASE; j++)
	  printf ("   Contador de %d: %ld\n", j, contador[j]);
	
	for (j = 0; j < BASE; j++)
	  contador[j] = 0;
	
	for (i = CESQUERDA + 1; i <= RESUMO + CESQUERDA; i++)
	  for (j = 0; j < BASE; j++)
	    if (pi.digitos[i] == j)
	      contador[j]++;
	
	printf ("\n");
	printf ("Estatistica: (%d casas, com %d casas inteiras e %d casas fracionarias\n", CASAS + CESQUERDA, 0, RESUMO);
	
	for (j = 0; j < BASE; j++)
	  printf ("   Contador de %d: %ld\n", j, contador[j]);
	
	
	return (0);
	}