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

#include "mpi.h"

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

#define DIM 10

int main(int argc, char *argv[])
{
	MPFVector x, b;
	MPFVector local_x, local_b;

	long int d_dim[10];
	mpf_t ip;
	long int i, j, local_dim;
	int myrank, num_procs;

	MPI_Init(&argc, &argv);
	MPI_Comm_rank(MPI_COMM_WORLD, &myrank);
	MPI_Comm_size(MPI_COMM_WORLD, &num_procs);

	_mpi_set_bnc_default_prec_decimal(50, MPI_COMM_WORLD);
	commit_mpf(&(MPI_MPF), ceil(50/log10(2.0)), MPI_COMM_WORLD);
	create_mpf_op(&(MPI_MPF_SUM), _mpi_mpf_add, MPI_COMM_WORLD);

	mpf_init(ip);

	local_dim = _mpi_divide_dim(d_dim, DIM, num_procs);
	if(myrank == 0)
	{
		x = init_mpfvector(local_dim * num_procs);
		b = init_mpfvector(local_dim * num_procs);

		for(i = 0; i < DIM; i++)
		{
			set_mpfvector_i_d(x, i, (double)(i + 1));
			set_mpfvector_i_d(b, i, (double)(DIM - i));
		}
	}
	local_x = _mpi_init_mpfvector(d_dim, DIM, MPI_COMM_WORLD);
	local_b = _mpi_init_mpfvector(d_dim, DIM, MPI_COMM_WORLD);

	_mpi_divide_mpfvector(local_b, d_dim, b, MPI_COMM_WORLD);
	_mpi_divide_mpfvector(local_x, d_dim, x, MPI_COMM_WORLD);

	_mpi_ip_mpfvector(ip, local_x, local_b, MPI_COMM_WORLD);

	_mpi_free_mpfvector(local_x);
	_mpi_free_mpfvector(local_b);

	/* free */
	if(myrank == 0)
	{
		printf("Inner Product(MPI):");
		mpf_out_str(stdout, 10, 0, ip);
		printf("\n");

		free_mpfvector(x);
		free_mpfvector(b);
	}

	free_mpf_op(&(MPI_MPF_SUM));
	free_mpf(&(MPI_MPF));
	MPI_Finalize();

	return EXIT_SUCCESS;
}