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

#include "mpi.h"

#define USE_GMP
#define USE_MPFR
#include "mpi_bnc.h"

#define MAX_LENGTH 1024

#define DEG 20

int main(int argc, char *argv[])
{
	long int i, dtimes, mpftimes;
	FILE * dcoef, * mpfcoef;

	long int local_dim, dd_dim[MPI_GMP_MAXPROCS];
	double startwtime[2], endwtime[2];
	CMPFArray cmpfans, cmpfinit, local_cmpfans, local_cmpfinit;
	MPFPoly mpff;
	mpf_t mpfabs_eps, mpfrel_eps;
	int myrank, num_procs;

/* for MPI */
	MPI_Init(&argc, &argv);
	MPI_Comm_rank(MPI_COMM_WORLD, &myrank);
	MPI_Comm_size(MPI_COMM_WORLD, &num_procs);

#define PREC 128
//#define PREC 256
//#define PREC 512
//#define PREC 1024
#define MPFDEG 20
//	set_bnc_default_prec(PREC);
	_mpi_set_bnc_default_prec(PREC, MPI_COMM_WORLD);
	commit_mpf(&(MPI_MPF), PREC, MPI_COMM_WORLD);
	commit_mpi_mpfcmplx(&(MPI_BNC_MPFCMPLX), PREC, MPI_COMM_WORLD);

	/* init */
	mpf_init_set_d(mpfabs_eps, 1.0e-100);
	mpf_init_set_d(mpfrel_eps, 1.0e-7);

	mpff = init_mpfpoly(MAX_LENGTH);
	cmpfans = init_cmpfarray(MPFDEG);
	cmpfinit = init_cmpfarray(MPFDEG);

	local_dim = _mpi_divide_dim(dd_dim, MPFDEG, num_procs);
	if(myrank == 0)
	{
//		for(i = 0; i < num_procs; i++)
//			printf("%d d_dim[%d]: %d\n", local_dim, i, dd_dim[i]);
	}
	local_cmpfans = init_cmpfarray(local_dim);
	local_cmpfinit = init_cmpfarray(local_dim);

	cmpfans = init_cmpfarray(MPFDEG);
	cmpfinit = init_cmpfarray(MPFDEG);
	if(myrank == 0)
	{
		mpfcoef = fopen("polycoef20.dat", "r");
		fread_mpfpolycoef(mpfcoef, mpff, DEG);
		print_mpfpoly(mpff);fflush(stdout);
		fclose(mpfcoef);
	}
//	goto end;

	/* Bcast mpfpoly */
	_mpi_bcast_mpfpoly(mpff, MPI_COMM_WORLD);
	if(myrank == 0)
	{
		printf("rank: %d\n", myrank);
		print_mpfpoly(mpff);fflush(stdout);
	}

	/* set Aberth's initial value */
	_mpi_mpf_dka_init(local_cmpfinit, mpff, MPI_COMM_WORLD);
//	print_cmpfarray(local_cmpfinit);

	/* DKA method */
	if(myrank == 0) startwtime[1] = MPI_Wtime();
	_mpi_mpf_dka(&mpftimes, cmpfans, local_cmpfans, cmpfinit, local_cmpfinit, mpff, 1000, mpfabs_eps, mpfrel_eps, MPI_COMM_WORLD);
	if(myrank == 0) endwtime[1] = MPI_Wtime();

	/* print answer */
	if(myrank == 0) printf("Iterative times: %d\n", mpftimes);
	_mpi_collect_cmpfarray(cmpfans, dd_dim, local_cmpfans, MPI_COMM_WORLD);
//	print_cmpfarray(local_cmpfans);
	if(myrank == 0) print_cmpfarray(cmpfans);

	/* clear */
	free_mpfpoly(mpff);
	free_cmpfarray(local_cmpfans);
	free_cmpfarray(local_cmpfinit);
	free_cmpfarray(cmpfans);
	free_cmpfarray(cmpfinit);

	free_mpi_mpfcmplx(&(MPI_BNC_MPFCMPLX));
	free_mpf(&(MPI_MPF));

end:
	MPI_Finalize();
	if(myrank == 0)
	{
		printf("mpf_DKA(%dbits): %f sec\n", PREC, endwtime[1] - startwtime[1]);
	}

	return EXIT_SUCCESS;

}