git » metnum-tp1.git » master » tree

[master] / binomial.cpp

#include <stdio.h>
#include <stdlib.h>
#include <stdint.h>
#include "aritmetica.hpp"
#include "util.hpp"

int prec;

static numero factorial(unsigned int n)
{
	numero res(prec,1);
	for (unsigned int i = 1; i <= n; i++) {
		res = res * i;
	}
	return res;
}

int main(int argc, char *argv[])
{
	if (argc != 3) {
		printf("Uso: %s <prec> <nterm>\n", argv[0]);
		return 1;
	}

	prec = atoi(argv[1]);
	int nterm = atoi(argv[2]);

	/* Valor inicial, el del primer término */
	numero n(prec, 1);
	print_result(1, n);

	/* Loop de la serie, nterm - 1 términos (el primero, por el cual
	 * arrancamos, ya esta fijo arriba) */
	for (int i = 1; i < nterm; i++) {

		/* numerador: 1/2 (1/2 - 1) (1/2 - 2) ...
		 * hasta (1/2 - i - 1) */
		numero unmedio(prec, 0.5);
		numero numerador = unmedio;
		for (double j = 1; j < i; j++) {
			numerador = numerador * (unmedio - j);
		}

		/* denominador: factorial(i) */
		numero denominador = factorial(i);

		/* Acumulamos el término */
		n = n + numerador / denominador;

		/* XXX DEBUG - SACAR */
		#if 0
		printf("  %.60f 0x%" print_64t "x\n",
				numerador.to_double(), to_u64(numerador));
		printf("  %.60f 0x%" print_64t "x\n",
				denominador.to_double(), to_u64(denominador));
		printf("  ---\n");
		#endif
		//printf(" %d %.60f 0x%" print_64t "x\n", i, n.to_double(), to_u64(n));

		/* Vamos emitiendo los resultados a medida que los tenemos */
		print_result(i + 1, n);
	}

	return 0;
}